source("/Shared_win/projects/RNA_normal/analysis.10x.r")
GEX.seur <- readRDS("./sn10x_SZJ.sct_anno.rds")
GEX.seur
## An object of class Seurat
## 47749 features across 28452 samples within 3 assays
## Active assay: SCT (21008 features, 0 variable features)
## 2 other assays present: RNA, integrated
## 3 dimensional reductions calculated: pca, tsne, umap
DefaultAssay(GEX.seur) <- 'integrated'
DefaultAssay(GEX.seur)
## [1] "integrated"
GEX.seur <- subset(GEX.seur, subset= seurat_clusters %in% c(0:29))
GEX.seur <- subset(GEX.seur, subset= percent.mt <1.5)
GEX.seur <- subset(GEX.seur, subset= (cnt %in% c("Stst.CTL","Stst.CKO") & nFeature_RNA > 500 & nFeature_RNA < 1800 & nCount_RNA < 4000) |
(cnt %in% c("CR7d.CTL","CR7d.CKO") & nFeature_RNA > 800 & nFeature_RNA < 2100 & nCount_RNA < 5000))
GEX.seur
## An object of class Seurat
## 47749 features across 24358 samples within 3 assays
## Active assay: integrated (2386 features, 2386 variable features)
## 2 other assays present: RNA, SCT
## 3 dimensional reductions calculated: pca, tsne, umap
length(GEX.seur@assays$integrated@var.features)
## [1] 2386
GEX.seur <- ScaleData(object = GEX.seur, verbose = TRUE,
vars.to.regress = c("percent.mt","percent.rb","nCount_RNA"))
## Regressing out percent.mt, percent.rb, nCount_RNA
## Centering and scaling data matrix
GEX.seur <- RunPCA(GEX.seur, do.print = TRUE,
features = GEX.seur@assays$integrated@var.features,
seed.use = 133,
npcs = 100,
#ndims.print = 1,
verbose = T)
## PC_ 1
## Positive: Nrg3, Grid2, Lrrtm4, Pcdh15, Trpm3, Nrg3os, Ndst4, Cadm1, Kcnq3, Nkain2
## Galntl6, Cacnb2, Slit2, Ctnna3, Frmpd4, Epha5, Magi1, Syn3, Kcnma1, Fstl5
## Pde4b, Kcnc2, Plce1, Grip1, Tox, Kcnd2, Ryr2, Sntb1, Grik2, Lrp1b
## Negative: Ntrk3, Robo2, Tmeff2, Nrxn3, Ano2, Cdh8, Cpne4, Mgat4c, Myl1, Clstn2
## Plxna4, Pdzrn4, Zfp804a, Adgrg6, Spock3, Dgkg, Gpr149, Cdh6, Pcdh10, Cux2
## Astn2, Cntn5, Cysltr2, Kcnb2, Iqgap2, Arhgap6, Grin3a, Cacna2d3, Ccbe1, Itgb6
## PC_ 2
## Positive: Rbfox1, Bnc2, Ptprt, Gpc6, Tafa1, Tshz2, Grik1, Mdga2, St6galnac3, Tox
## Adgrb3, Frmd4b, Brinp2, Fbxw15, Pcdh7, Cdc14a, Plcxd3, Agbl4, Fbxw24, Pbx1
## Unc5c, Chat, Pld5, Oprk1, Negr1, Dock2, Pde4b, Gfra2, Adamtsl1, Caln1
## Negative: Nos1, Auts2, Etv1, Cadps2, Gfra1, Alcam, Fam155a, Kcnab1, Egfem1, Kcnq5
## Asic2, Dgkb, Dach1, Plekha5, Schip1, Rgs6, Kcnt2, Epha5, Cmah, Ank2
## Nav3, Stxbp6, Cntnap5a, Creb5, Hs6st3, Tmem108, Lrrc4c, Ncam2, Frmpd4, Ablim2
## PC_ 3
## Positive: Kcnip4, Cdh18, Csmd3, Kctd8, Cadm2, Klhl1, Gabrg3, Cntn3, Pbx3, Htr4
## Pde4d, Dlc1, Prkg1, Dmd, Dock10, Khdrbs2, Car10, Edil3, Meis1, Gpc6
## Skap1, Tenm2, Serpini1, Tac1, Plcl1, Grm7, C79798, Gda, March1, Nrp2
## Negative: Sgcd, Ptprg, Adgrg6, Nfia, Filip1, Cysltr2, Ccbe1, Slc4a4, Nfib, Gpr149
## Fgf13, Nmu, Sema6d, Ano2, Dgkg, Grin3a, Cpne4, Dapk2, Itgb6, Malat1
## Zfp804a, Ngfr, Cbln2, Efr3a, Nos1, Tshz2, Lhfpl2, Airn, Cntnap5a, Zfp521
## PC_ 4
## Positive: Kcnt2, Lingo2, Fgf13, Prkg1, Dock10, Ndst4, Epha5, Lrrc4c, Ctnna3, Cntn5
## Lrrtm4, Dmd, Gda, Kcnip4, Tac1, Chl1, Ank2, Sorcs1, Rgs6, Nxph1
## Ptprz1, Hs3st2, Thsd4, Necab1, Grem2, Rora, Kctd8, Man1a, Galntl6, Sorcs3
## Negative: Trhde, Chsy3, Ebf1, Trps1, March1, Gal, Cntn4, Col18a1, Nrg1, Enox1
## Trpm3, Zmat4, Ntng1, Cpa6, Sdk1, Csmd1, Dcc, Shisa6, Ccser1, Tenm4
## Npas3, Nkain3, Kcnd2, Plcxd3, Tenm1, Sctr, Nwd2, Kcnh7, Prune2, Sez6l
## PC_ 5
## Positive: Ptprd, Nrg1, Trhde, Lsamp, Egfem1, Cntn4, Cntn3, Adgrl2, Rmst, Cntn5
## Csmd3, Sgcz, Kcnd2, Luzp2, Gal, Nav2, Astn2, Ebf1, Asic2, Trps1
## Cpa6, Lrp1b, Kcnip4, Hs6st3, Sorcs1, Zbtb7c, Moxd1, Scn11a, Zmat4, Csmd1
## Negative: Dgkb, Klhl1, Vwc2l, Pbx3, Il1rapl1, Rasgef1b, Cdh9, Alk, Zfhx3, Mgat4c
## Sema5a, Dpp10, Vcan, Alcam, Galnt18, Galnt13, Auts2, Fam155a, Zbbx, C79798
## Olfr78, Thsd7b, Pcdh7, Scgn, Nek1, Serpini1, Sntg1, Lncbate1, P3h2, Mir100hg
DimHeatmap(GEX.seur, dims = 1:12, cells = 1500, balanced = TRUE,ncol = 4)
ElbowPlot(GEX.seur, ndims = 100)
ElbowPlot(GEX.seur, ndims = 50)
#
PCsct <- 1:18
GEX.seur <- FindNeighbors(GEX.seur, k.param = 20, dims = PCsct, compute.SNN = T, reduction = 'pca', verbose = T)
## Computing nearest neighbor graph
## Computing SNN
GEX.seur <- FindClusters(GEX.seur, dims.use = PCsct, algorithm = 1, save.SNN =T, resolution = 2, reduction = 'pca', verbose = T)
## Warning: The following arguments are not used: dims.use, save.SNN, reduction
## Suggested parameter: dims instead of dims.use
## Warning: The following arguments are not used: dims.use, save.SNN, reduction
## Suggested parameter: dims instead of dims.use
## Modularity Optimizer version 1.3.0 by Ludo Waltman and Nees Jan van Eck
##
## Number of nodes: 24358
## Number of edges: 917281
##
## Running Louvain algorithm...
## Maximum modularity in 10 random starts: 0.8329
## Number of communities: 26
## Elapsed time: 4 seconds
GEX.seur <- RunTSNE(object = GEX.seur, assay = "integrated", seed.use = 233, dims = PCsct, complexity=100)
GEX.seur <- RunUMAP(object = GEX.seur, assay = "integrated", seed.use = 166, dims = PCsct, n.neighbors = 20, min.dist = 0.3)
## Warning: The default method for RunUMAP has changed from calling Python UMAP via reticulate to the R-native UWOT using the cosine metric
## To use Python UMAP via reticulate, set umap.method to 'umap-learn' and metric to 'correlation'
## This message will be shown once per session
## 00:44:06 UMAP embedding parameters a = 0.9922 b = 1.112
## 00:44:06 Read 24358 rows and found 18 numeric columns
## 00:44:06 Using Annoy for neighbor search, n_neighbors = 20
## 00:44:06 Building Annoy index with metric = cosine, n_trees = 50
## 0% 10 20 30 40 50 60 70 80 90 100%
## [----|----|----|----|----|----|----|----|----|----|
## **************************************************|
## 00:44:08 Writing NN index file to temp file C:\Users\Shaorui\AppData\Local\Temp\RtmpMrRNx9\file848f084ce9
## 00:44:08 Searching Annoy index using 1 thread, search_k = 2000
## 00:44:14 Annoy recall = 100%
## 00:44:14 Commencing smooth kNN distance calibration using 1 thread with target n_neighbors = 20
## 00:44:16 Initializing from normalized Laplacian + noise (using irlba)
## 00:44:18 Commencing optimization for 200 epochs, with 703756 positive edges
## 00:44:42 Optimization finished
DimPlot(GEX.seur, label = T, pt.size = 0.05, repel = F, reduction = 'tsne', group.by = "seurat_clusters") +
DimPlot(GEX.seur, label = T, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "seurat_clusters")
DimPlot(GEX.seur, label = T, pt.size = 0.05, repel = F, reduction = 'tsne', group.by = "Anno1") +
DimPlot(GEX.seur, label = T, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "Anno1")
DimPlot(GEX.seur, label = T, pt.size = 0.05, repel = F, reduction = 'tsne', group.by = "Anno2") +
DimPlot(GEX.seur, label = T, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "Anno2")
color.cnt <- scales::hue_pal()(4)[c(2,1,3,4)]
color.cnt
## [1] "#7CAE00" "#F8766D" "#00BFC4" "#C77CFF"
DimPlot(GEX.seur, label = F, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "cnt", split.by = "cnt", ncol = 4, cols = color.cnt)
DimPlot(GEX.seur, label = T, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "Anno1") +
DimPlot(GEX.seur, label = F, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "cnt",cols = color.cnt)
DimPlot(subset(GEX.seur, subset = cnt %in% c("Stst.CTL","Stst.CKO")), label = F, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "cnt", cols = color.cnt[1:2]) +
DimPlot(subset(GEX.seur, subset = cnt %in% c("CR7d.CTL","CR7d.CKO")), label = F, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "cnt", cols = color.cnt[3:4])
FeaturePlot(GEX.seur,
reduction = "umap", features = c("nFeature_RNA","nCount_RNA","percent.mt","percent.rb"))
FeaturePlot(subset(GEX.seur, subset = cnt %in% c("Stst.CTL","Stst.CKO")),
reduction = "umap", features = c("nFeature_RNA","nCount_RNA","percent.mt","percent.rb"))
FeaturePlot(subset(GEX.seur, subset = cnt %in% c("CR7d.CTL","CR7d.CKO")),
reduction = "umap", features = c("nFeature_RNA","nCount_RNA","percent.mt","percent.rb"))
(DimPlot(subset(GEX.seur,subset=cnt %in% c("Stst.CTL","Stst.CKO")), label = T, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "Anno1")+labs(title="Stst only")) +
(DimPlot(subset(GEX.seur,subset=cnt %in% c("CR7d.CTL","CR7d.CKO")), label = T, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "Anno1")+labs(title="CR7d only"))
(DimPlot(subset(GEX.seur,subset=cnt %in% c("Stst.CTL","Stst.CKO")), label = T, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "seurat_clusters")+labs(title="Stst only")) +
(DimPlot(subset(GEX.seur,subset=cnt %in% c("CR7d.CTL","CR7d.CKO")), label = T, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "seurat_clusters")+labs(title="CR7d only"))
(DimPlot(subset(GEX.seur,subset=cnt %in% c("Stst.CTL","Stst.CKO")), label = T, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "Anno1")+labs(title="Stst only")) +
(DimPlot(subset(GEX.seur,subset=cnt %in% c("Stst.CTL","Stst.CKO")), label = T, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "seurat_clusters")+labs(title="Stst only"))
(DimPlot(subset(GEX.seur,subset=cnt %in% c("CR7d.CTL","CR7d.CKO")), label = T, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "Anno1")+labs(title="CR7d only")) +
(DimPlot(subset(GEX.seur,subset=cnt %in% c("CR7d.CTL","CR7d.CKO")), label = T, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "seurat_clusters")+labs(title="CR7d only"))
GEX.seur$sort_clusters <- factor(as.character(GEX.seur$seurat_clusters),
levels = c(3,0,7,8, 14, 13, 9, 1, 17,
6,5,2,4, 12, 16, 23,
10, 19,25, 22,
18,11, 15,21, 24, 20))
VlnPlot(subset(GEX.seur, subset = cnt %in% c("Stst.CTL","Stst.CKO")),
features = c("nFeature_RNA", "nCount_RNA", "percent.mt","percent.rb"), ncol = 2, pt.size = 0.01, group.by = "sort_clusters")
VlnPlot(subset(GEX.seur, subset = cnt %in% c("CR7d.CTL","CR7d.CKO")),
features = c("nFeature_RNA", "nCount_RNA", "percent.mt","percent.rb"), ncol = 2, pt.size = 0.01, group.by = "sort_clusters")
DefaultAssay(GEX.seur) <- "SCT"
DefaultAssay(GEX.seur)
## [1] "SCT"
markers.new.ss <- list(EMN=c("Chat","Bnc2","Tox","Ptprt",
"Gfra2","Oprk1","Adamtsl1",
"Fbxw15","Fbxw24","Chrna7",
"Satb1","Itga6","Cntnap5b",
"Pgm5","Chgb",
"Nxph1",
"Gabrb1","Glp2r","Nebl",
"Lrrc7",
"Ryr3","Eda",
"Hgf","Lama2","Efnb2",
"Tac1",
"Kctd8",
"Ptn",
"Ntrk2","Penk","Itgb8",
"Fut9","Nfatc1","Egfr","Pdia5",
"Ahr","Mgll",
"Aff3",
"Chrm3"
),
IMN=c("Nos1","Kcnab1",
"Gfra1","Etv1",
"Man1a","Airn",
"Adcy2",
"Col25a1",
"Cmah","Creb5","Vip","Pde1a",
"Ebf1","Gpc5","Mid1","Igfbp5",
"Ppara",
"Pcdh11x","Adcy8","Grp"
),
IN=c("Npas3","Synpr","St18","Gal",
"Nova1",
"Cdh10","Kcnk13",
"Neurod6","Moxd1","Sctr",
"Piezo1","Vipr2","Adamts9","Sst","Kcnn2"
),
IPAN=c("Calb2","Adcy1","Calcb","Nmu","Adgrg6","Pcdh10",
"Ngfr","Galr1","Il7","Aff2",
"Gpr149","Itgb6","Met","Itgbl1",
"Cdh6","Cdh8",
"Clstn2","Ano2","Ntrk3",
"Cpne4","Vwc2l","Cdh9",
"Car10","Dcc",
"Scgn",
"Vcan","Cck","Piezo2","Kcnh7",
"Rerg","Bmpr1b","Skap1","Ntng1",
"Tafa2","Nxph2"),
CONTAM=c("Actl6b","Phox2b"), #,"Sox10","Plp1","L1cam","Gfap","Rxrg","Acta2","Actg2","Epcam","Pecam1","Ptprc"
pDEGs=c("Grid2","Ide","Btaf1","Slc15a2","Ccser1","Hdac9","Rspo2","Grem2"))
#
pm.All_new.s2 <- DotPlot(GEX.seur, features = as.vector(unlist(markers.new.ss)), group.by = "sort_clusters",
cols = c("midnightblue","darkorange1")) +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+ scale_y_discrete(limits=rev) #+ labs(title="All - newAnno1")
pm.All_new.s2
C25 ?
table(GEX.seur@meta.data[,c("sample","sort_clusters")])
## sort_clusters
## sample 3 0 7 8 14 13 9 1 17 6 5 2 4 12 16 23 10
## CR7d.CKO1 147 66 176 100 60 72 70 123 41 116 88 102 112 45 47 16 111
## CR7d.CKO2 86 38 111 68 51 41 64 83 26 69 46 66 91 43 32 14 73
## CR7d.CKO3 107 38 155 113 67 76 73 100 54 105 64 113 116 52 61 22 93
## CR7d.CTL1 127 72 223 133 74 74 97 134 76 145 90 132 122 63 57 20 71
## CR7d.CTL2 141 60 199 136 93 90 96 143 52 141 94 115 119 71 45 21 57
## CR7d.CTL3 118 43 192 105 58 50 71 113 49 121 90 126 126 82 59 28 98
## CR7d.CTL4 150 63 188 149 74 72 104 141 48 155 83 126 119 95 60 20 85
## Stst.CKO1 95 244 4 45 69 68 93 135 53 58 132 146 92 88 48 44 129
## Stst.CKO2 137 305 4 55 75 102 96 137 60 61 123 156 108 71 63 33 106
## Stst.CKO3 99 224 2 64 66 79 110 99 54 62 113 110 76 62 61 22 93
## Stst.CTL1 135 277 8 84 68 72 104 169 55 100 139 152 74 115 54 26 78
## Stst.CTL2 93 174 6 71 46 55 100 111 59 68 106 112 89 102 42 24 69
## Stst.CTL3 128 306 10 62 72 70 81 144 53 78 134 136 79 88 53 22 63
## sort_clusters
## sample 19 25 22 18 11 15 21 24 20
## CR7d.CKO1 42 2 22 47 75 41 45 8 49
## CR7d.CKO2 30 6 23 34 70 40 29 9 20
## CR7d.CKO3 42 9 37 31 75 52 37 11 40
## CR7d.CTL1 48 4 19 41 101 62 41 15 42
## CR7d.CTL2 29 5 31 23 72 49 51 21 34
## CR7d.CTL3 60 5 28 40 67 45 60 10 36
## CR7d.CTL4 46 4 42 48 97 67 48 22 48
## Stst.CKO1 55 3 30 63 77 66 30 18 40
## Stst.CKO2 51 5 17 65 98 80 25 20 51
## Stst.CKO3 39 1 21 47 64 71 34 14 46
## Stst.CTL1 51 1 15 44 65 66 35 15 42
## Stst.CTL2 39 2 16 39 83 65 25 11 38
## Stst.CTL3 39 3 16 67 90 70 35 17 57
I perfer to remove C25 which might be a rare new subtype, but here more like a contaminated mix in IN2
GEX.seur <- subset(GEX.seur, subset = seurat_clusters != c(25))
GEX.seur
## An object of class Seurat
## 47749 features across 24308 samples within 3 assays
## Active assay: SCT (21008 features, 0 variable features)
## 2 other assays present: RNA, integrated
## 3 dimensional reductions calculated: pca, tsne, umap
#
GEX.seur@meta.data[,grep("Anno|snn",colnames(GEX.seur@meta.data),value = T)] <- NULL
head(GEX.seur@meta.data)
## orig.ident nCount_RNA nFeature_RNA percent.mt percent.rb
## AAACCCAAGAATACAC-1_2 Stst.CTL_CKO 2282 1290 0.1314636 0.1314636
## AAACCCAAGCAATAGT-1_2 Stst.CTL_CKO 1535 1011 0.1302932 0.5863192
## AAACCCAAGGTGAGCT-1_2 Stst.CTL_CKO 1049 737 0.0000000 0.4766444
## AAACCCACAACGAGGT-1_2 Stst.CTL_CKO 1591 1046 0.6913891 0.4399749
## AAACCCACAAGAGTAT-1_2 Stst.CTL_CKO 2106 1250 0.2849003 0.2849003
## AAACCCACAATCAAGA-1_2 Stst.CTL_CKO 1451 953 0.1378360 0.1378360
## FB.info S.Score G2M.Score Phase cnt rep
## AAACCCAAGAATACAC-1_2 CTL.3 -0.00110357 0.012058898 G2M Stst.CTL rep3
## AAACCCAAGCAATAGT-1_2 CKO.3 0.02907907 -0.015558699 S Stst.CKO rep3
## AAACCCAAGGTGAGCT-1_2 CTL.1 -0.01100811 -0.008015087 G1 Stst.CTL rep1
## AAACCCACAACGAGGT-1_2 CTL.1 -0.02143685 0.005086498 G2M Stst.CTL rep1
## AAACCCACAAGAGTAT-1_2 CKO.3 0.02965845 -0.009057774 S Stst.CKO rep3
## AAACCCACAATCAAGA-1_2 CTL.2 -0.01216686 0.004143546 G2M Stst.CTL rep2
## sample tissue nCount_SCT nFeature_SCT seurat_clusters
## AAACCCAAGAATACAC-1_2 Stst.CTL3 Ileum 1860 1287 8
## AAACCCAAGCAATAGT-1_2 Stst.CKO3 Ileum 1525 1010 2
## AAACCCAAGGTGAGCT-1_2 Stst.CTL1 Ileum 1168 735 0
## AAACCCACAACGAGGT-1_2 Stst.CTL1 Ileum 1568 1045 3
## AAACCCACAAGAGTAT-1_2 Stst.CKO3 Ileum 1831 1248 16
## AAACCCACAATCAAGA-1_2 Stst.CTL2 Ileum 1451 953 3
## sort_clusters
## AAACCCAAGAATACAC-1_2 8
## AAACCCAAGCAATAGT-1_2 2
## AAACCCAAGGTGAGCT-1_2 0
## AAACCCACAACGAGGT-1_2 3
## AAACCCACAAGAGTAT-1_2 16
## AAACCCACAATCAAGA-1_2 3
# intAnno1
GEX.seur$intAnno1 <- as.character(GEX.seur$seurat_clusters)
GEX.seur$intAnno1[GEX.seur$intAnno1 %in% c(3,0,7,8,14)] <- "EMN1"
GEX.seur$intAnno1[GEX.seur$intAnno1 %in% c(13)] <- "EMN2"
GEX.seur$intAnno1[GEX.seur$intAnno1 %in% c(9)] <- "EMN3"
GEX.seur$intAnno1[GEX.seur$intAnno1 %in% c(1)] <- "EMN4"
GEX.seur$intAnno1[GEX.seur$intAnno1 %in% c(17)] <- "EMN5"
GEX.seur$intAnno1[GEX.seur$intAnno1 %in% c(6,5,2,4)] <- "IMN1"
GEX.seur$intAnno1[GEX.seur$intAnno1 %in% c(12)] <- "IMN2"
GEX.seur$intAnno1[GEX.seur$intAnno1 %in% c(16)] <- "IMN3"
GEX.seur$intAnno1[GEX.seur$intAnno1 %in% c(23)] <- "IMN4"
GEX.seur$intAnno1[GEX.seur$intAnno1 %in% c(10)] <- "IN1"
GEX.seur$intAnno1[GEX.seur$intAnno1 %in% c(19)] <- "IN2"
GEX.seur$intAnno1[GEX.seur$intAnno1 %in% c(22)] <- "IN3"
GEX.seur$intAnno1[GEX.seur$intAnno1 %in% c(18,11)] <- "IPAN1"
GEX.seur$intAnno1[GEX.seur$intAnno1 %in% c(15,21)] <- "IPAN2"
GEX.seur$intAnno1[GEX.seur$intAnno1 %in% c(24)] <- "IPAN3"
GEX.seur$intAnno1[GEX.seur$intAnno1 %in% c(20)] <- "IPAN4"
GEX.seur$intAnno1 <- factor(GEX.seur$intAnno1,
levels = c(paste0("EMN",1:5),
paste0("IMN",1:4),
paste0("IN",1:3),
paste0("IPAN",1:4)))
# intAnno2
GEX.seur$intAnno2 <- as.character(GEX.seur$seurat_clusters)
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(3,0,7,8,14)] <- "EMN1"
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(13)] <- "EMN2"
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(9)] <- "EMN3"
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(1)] <- "EMN4"
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(17)] <- "EMN5"
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(6,5,2,4)] <- "IMN1"
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(12)] <- "IMN2"
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(16)] <- "IMN3"
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(23)] <- "IMN4"
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(10)] <- "IN1"
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(19)] <- "IN2"
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(22)] <- "IN3"
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(18)] <- "IPAN1.1"
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(11)] <- "IPAN1.2"
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(15)] <- "IPAN2.1"
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(21)] <- "IPAN2.2"
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(24)] <- "IPAN3"
GEX.seur$intAnno2[GEX.seur$intAnno2 %in% c(20)] <- "IPAN4"
GEX.seur$intAnno2 <- factor(GEX.seur$intAnno2,
levels = c(paste0("EMN",1:5),
paste0("IMN",1:4),
paste0("IN",1:3),
paste0("IPAN",c(1.1,1.2,2.1,2.2,3,4))))
#saveRDS(GEX.seur, "./sn10x_SZJ.sct_anno.s.rds")
#GEX.seur <- readRDS("./sn10x_SZJ.sct_anno.s.rds")
pp.umap1 <- DimPlot(GEX.seur, label = T, pt.size = 0.05, repel = F, reduction = 'tsne', group.by = "seurat_clusters") +
DimPlot(GEX.seur, label = T, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "seurat_clusters")
pp.umap1
pp.umap2 <- DimPlot(GEX.seur, label = T, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "intAnno1") +
DimPlot(GEX.seur, label = T, label.size = 2.65, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "intAnno2")
pp.umap2
pp.umap3 <- DimPlot(GEX.seur, label = F, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "cnt", split.by = "cnt", ncol = 4, cols = color.cnt)
pp.umap3
markers.new.ss <- list(EMN=c("Chat","Bnc2","Tox","Ptprt",
"Gfra2","Oprk1","Adamtsl1",
"Fbxw15","Fbxw24","Chrna7",
"Satb1","Itga6","Cntnap5b",
"Pgm5","Chgb",
"Nxph1",
"Gabrb1","Glp2r","Nebl",
"Lrrc7",
"Ryr3","Eda",
"Hgf","Lama2","Efnb2",
"Tac1",
"Kctd8",
"Ptn",
"Ntrk2","Penk","Itgb8",
"Fut9","Nfatc1","Egfr","Pdia5",
"Ahr","Mgll",
"Aff3",
"Chrm3"
),
IMN=c("Nos1","Kcnab1",
"Gfra1","Etv1",
"Man1a","Airn",
"Adcy2",
"Col25a1",
"Cmah","Creb5","Vip","Pde1a",
"Ebf1","Gpc5","Mid1","Igfbp5",
"Ppara",
"Pcdh11x","Adcy8","Grp"
),
IN=c("Npas3","Synpr","St18","Gal",
"Nova1",
"Cdh10","Kcnk13",
"Neurod6","Moxd1","Sctr",
"Piezo1","Vipr2","Adamts9","Sst","Kcnn2"
),
IPAN=c("Calb2","Adcy1","Calcb","Nmu","Adgrg6","Pcdh10",
"Ngfr","Galr1","Il7","Aff2",
"Gpr149","Itgb6","Met","Itgbl1",
"Cdh6","Cdh8",
"Clstn2","Ano2","Ntrk3",
"Cpne4","Vwc2l","Cdh9",
"Car10","Dcc",
"Scgn",
"Vcan","Cck","Piezo2","Kcnh7",
"Rerg","Bmpr1b","Skap1","Ntng1",
"Tafa2","Nxph2"),
CONTAM=c("Actl6b","Phox2b"),
pDEGs=c("Grid2","Ide","Btaf1","Slc15a2","Ccser1","Hdac9","Rspo2","Grem2"))
check.markers.ss <- as.vector(unlist(markers.new.ss))
check.markers.ss
## [1] "Chat" "Bnc2" "Tox" "Ptprt" "Gfra2" "Oprk1"
## [7] "Adamtsl1" "Fbxw15" "Fbxw24" "Chrna7" "Satb1" "Itga6"
## [13] "Cntnap5b" "Pgm5" "Chgb" "Nxph1" "Gabrb1" "Glp2r"
## [19] "Nebl" "Lrrc7" "Ryr3" "Eda" "Hgf" "Lama2"
## [25] "Efnb2" "Tac1" "Kctd8" "Ptn" "Ntrk2" "Penk"
## [31] "Itgb8" "Fut9" "Nfatc1" "Egfr" "Pdia5" "Ahr"
## [37] "Mgll" "Aff3" "Chrm3" "Nos1" "Kcnab1" "Gfra1"
## [43] "Etv1" "Man1a" "Airn" "Adcy2" "Col25a1" "Cmah"
## [49] "Creb5" "Vip" "Pde1a" "Ebf1" "Gpc5" "Mid1"
## [55] "Igfbp5" "Ppara" "Pcdh11x" "Adcy8" "Grp" "Npas3"
## [61] "Synpr" "St18" "Gal" "Nova1" "Cdh10" "Kcnk13"
## [67] "Neurod6" "Moxd1" "Sctr" "Piezo1" "Vipr2" "Adamts9"
## [73] "Sst" "Kcnn2" "Calb2" "Adcy1" "Calcb" "Nmu"
## [79] "Adgrg6" "Pcdh10" "Ngfr" "Galr1" "Il7" "Aff2"
## [85] "Gpr149" "Itgb6" "Met" "Itgbl1" "Cdh6" "Cdh8"
## [91] "Clstn2" "Ano2" "Ntrk3" "Cpne4" "Vwc2l" "Cdh9"
## [97] "Car10" "Dcc" "Scgn" "Vcan" "Cck" "Piezo2"
## [103] "Kcnh7" "Rerg" "Bmpr1b" "Skap1" "Ntng1" "Tafa2"
## [109] "Nxph2" "Actl6b" "Phox2b" "Grid2" "Ide" "Btaf1"
## [115] "Slc15a2" "Ccser1" "Hdac9" "Rspo2" "Grem2"
length(check.markers.ss)
## [1] 119
sum(check.markers.ss %in% rownames(GEX.seur))
## [1] 119
pm.All_new.a <- DotPlot(GEX.seur, features = check.markers.ss, group.by = "sort_clusters",
cols = c("midnightblue","darkorange1")) +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+ scale_y_discrete(limits=rev) + labs(title="All - sort_clusters")
pm.All_new.a
pm.All_new.b <- DotPlot(GEX.seur, features = check.markers.ss, group.by = "intAnno1",
cols = c("midnightblue","darkorange1")) +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+ scale_y_discrete(limits=rev) + labs(title="All - intAnno1")
pm.All_new.b
pm.All_new.c <- DotPlot(GEX.seur, features = check.markers.ss, group.by = "intAnno2",
cols = c("midnightblue","darkorange1")) +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+ scale_y_discrete(limits=rev) + labs(title="All - intAnno2")
pm.All_new.c
#
pm.All_new.c1 <- DotPlot(subset(GEX.seur, subset= cnt %in% c("Stst.CTL")), features = check.markers.ss, group.by = "intAnno2",
cols = c("midnightblue","darkorange1")) +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+ scale_y_discrete(limits=rev) + labs(title="Stst.CTL - intAnno2")
pm.All_new.c1
#
pm.All_new.c2 <- DotPlot(subset(GEX.seur, subset= cnt %in% c("Stst.CKO")), features = check.markers.ss, group.by = "intAnno2",
cols = c("midnightblue","darkorange1")) +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+ scale_y_discrete(limits=rev) + labs(title="Stst.CKO - intAnno2")
pm.All_new.c2
#
pm.All_new.c3 <- DotPlot(subset(GEX.seur, subset= cnt %in% c("CR7d.CTL")), features = check.markers.ss, group.by = "intAnno2",
cols = c("midnightblue","darkorange1")) +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+ scale_y_discrete(limits=rev) + labs(title="CR7d.CTL - intAnno2")
pm.All_new.c3
#
pm.All_new.c4 <- DotPlot(subset(GEX.seur, subset= cnt %in% c("CR7d.CKO")), features = check.markers.ss, group.by = "intAnno2",
cols = c("midnightblue","darkorange1")) +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+ scale_y_discrete(limits=rev) + labs(title="CR7d.CKO - intAnno2")
pm.All_new.c4
# define intAnno1/2 colors
color.A1 <- c("#678BB1","#8AB6CE","#3975C1","#669FDF","#4CC1BD",
"#BF7E6B","#D46B35","#F19258","#FF8080",
"#BDAE8D","#BD66C4","#C03778",
"#97BA59","#DFAB16","#2BA956","#9FE727")
names(color.A1) <- levels(GEX.seur$intAnno1)
color.A2 <- c("#678BB1","#8AB6CE","#3975C1","#669FDF","#4CC1BD",
"#BF7E6B","#D46B35","#F19258","#FF8080",
"#BDAE8D","#BD66C4","#C03778",
"#97BA59","#C4D116", "#DFAB16","#EDE25A", "#2BA956","#9FE727")
names(color.A2) <- levels(GEX.seur$intAnno2)
# define batch/condition colors
color.cnt <- scales::hue_pal()(4)[c(2,1,3,4)]
names(color.cnt) <- levels(GEX.seur$cnt)
color.test1 <- color.cnt[1:2]
color.test2 <- color.cnt[3:4]
## define feature colors
# Cell2020
material.heat = function(n){
colorRampPalette(
c(
"#283593", # indigo 800
"#3F51B5", # indigo
"#2196F3", # blue
"#00BCD4", # cyan
"#4CAF50", # green
"#8BC34A", # light green
"#CDDC39", # lime
"#FFEB3B", # yellow
"#FFC107", # amber
"#FF9800", # orange
"#FF5722" # deep orange
)
)(n)
}
# Immunity2019, na gray
colors.Immunity <-c("#191970","#121285","#0C0C9A","#0707B0","#0101C5","#0014CF","#0033D3","#0053D8","#0072DD","#0092E1","#00B2E6",
"#00D1EB","#23E8CD","#7AF17B","#D2FA29","#FFEB00","#FFC300","#FF9B00","#FF8400","#FF7800","#FF6B00","#FF5F00","#FF5300",
"#FF4700","#F73B00","#EF2E00","#E62300","#DD1700","#D50B00","#CD0000")
# NatNeur2021, sc-neurons
color.ref <- c("#8AB6CE","#678BB1","#3975C1","#4CC1BD",
"#C03778","#97BA59","#DFAB16","#BF7E6B",
"#D46B35","#BDAE8D","#BD66C4","#2BA956",
"#FF8080","#FF8080","#FF8080","#FF0000")
cowplot::plot_grid(
DimPlot(GEX.seur, reduction = "umap", group.by = "intAnno1", label = T, label.size = 3.25,repel = F, pt.size = 0.05,
cols = color.A1),
DimPlot(GEX.seur, label = F, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "cnt",
cols = color.cnt) ,
rel_widths = c(4.8,5),
ncol = 2)
cowplot::plot_grid(
DimPlot(GEX.seur, reduction = "umap", group.by = "intAnno2", label = T, label.size = 2.65,repel = F, pt.size = 0.05,
cols = color.A2),
DimPlot(GEX.seur, label = F, pt.size = 0.05, repel = F, reduction = 'umap', group.by = "cnt",
cols = color.cnt) ,
rel_widths = c(4.85,5),
ncol = 2)
test1.seur <- subset(GEX.seur, subset= cnt %in% c("Stst.CTL","Stst.CKO"))
test1.seur
## An object of class Seurat
## 47749 features across 11409 samples within 3 assays
## Active assay: SCT (21008 features, 0 variable features)
## 2 other assays present: RNA, integrated
## 3 dimensional reductions calculated: pca, tsne, umap
cowplot::plot_grid(
DimPlot(test1.seur, reduction = "umap", group.by = "intAnno1", label = T, label.size = 3.25,repel = F, pt.size = 0.15,
cols = color.A1),
DimPlot(test1.seur, label = F, pt.size = 0.15, repel = F, reduction = 'umap', group.by = "cnt",
cols = color.test1) ,
rel_widths = c(4.8,5),
ncol = 2)
cowplot::plot_grid(
DimPlot(test1.seur, reduction = "umap", group.by = "intAnno2", label = T, label.size = 2.65,repel = F, pt.size = 0.15,
cols = color.A2),
DimPlot(test1.seur, label = F, pt.size = 0.15, repel = F, reduction = 'umap', group.by = "cnt",
cols = color.test1) ,
rel_widths = c(4.85,5),
ncol = 2)
test2.seur <- subset(GEX.seur, subset= cnt %in% c("CR7d.CTL","CR7d.CKO"))
test2.seur
## An object of class Seurat
## 47749 features across 12899 samples within 3 assays
## Active assay: SCT (21008 features, 0 variable features)
## 2 other assays present: RNA, integrated
## 3 dimensional reductions calculated: pca, tsne, umap
cowplot::plot_grid(
DimPlot(test2.seur, reduction = "umap", group.by = "intAnno1", label = T, label.size = 3.25,repel = F, pt.size = 0.15,
cols = color.A1),
DimPlot(test2.seur, label = F, pt.size = 0.15, repel = F, reduction = 'umap', group.by = "cnt",
cols = color.test2) ,
rel_widths = c(4.8,5),
ncol = 2)
cowplot::plot_grid(
DimPlot(test2.seur, reduction = "umap", group.by = "intAnno2", label = T, label.size = 2.65,repel = F, pt.size = 0.15,
cols = color.A2),
DimPlot(test2.seur, label = F, pt.size = 0.15, repel = F, reduction = 'umap', group.by = "cnt",
cols = color.test2) ,
rel_widths = c(4.85,5),
ncol = 2)
cowplot::plot_grid(
pheatmap::pheatmap(table(cnt=test1.seur$sample,
clusters=test1.seur$intAnno1)[c(4:6,1:3),],
main = "Cell Count",
gaps_row = c(3),
gaps_col = c(5,9,12),
cluster_rows = F, cluster_cols = F, display_numbers = T, number_format = "%.0f", legend = F, silent = T)$gtable,
pheatmap::pheatmap(100*rowRatio(table(cnt=test1.seur$sample,
clusters=test1.seur$intAnno1)[c(4:6,1:3),] ),
main = "Cell Ratio",
gaps_row = c(3),
gaps_col = c(5,9,12),
cluster_rows = F, cluster_cols = F, display_numbers = T, number_format = "%.1f", legend = F, silent =T)$gtable,
ncol = 1)
cowplot::plot_grid(
pheatmap::pheatmap(table(cnt=test1.seur$sample,
clusters=test1.seur$intAnno2)[c(4:6,1:3),],
main = "Cell Count",
gaps_row = c(3),
gaps_col = c(5,9,12),
cluster_rows = F, cluster_cols = F, display_numbers = T, number_format = "%.0f", legend = F, silent = T)$gtable,
pheatmap::pheatmap(100*rowRatio(table(cnt=test1.seur$sample,
clusters=test1.seur$intAnno2)[c(4:6,1:3),] ),
main = "Cell Ratio",
gaps_row = c(3),
gaps_col = c(5,9,12),
cluster_rows = F, cluster_cols = F, display_numbers = T, number_format = "%.1f", legend = F, silent =T)$gtable,
ncol = 1)
stat1_intAnno1 <- test1.seur@meta.data[,c("cnt","rep","intAnno1")]
stat1_intAnno1.s <- stat1_intAnno1 %>%
group_by(cnt,rep,intAnno1) %>%
summarise(count=n()) %>%
mutate(prop= count/sum(count)) %>%
ungroup
## `summarise()` has grouped output by 'cnt', 'rep'. You can override using the
## `.groups` argument.
pcom1.intAnno1 <- stat1_intAnno1.s %>%
ggplot(aes(x = intAnno1, y = 100*prop, color=cnt)) +
geom_bar(stat="summary", fun="mean", position = position_dodge(0.75), width = 0.58, fill="white") +
theme_classic(base_size = 15) +
scale_color_manual(values = color.test1, name="") +
labs(title="Cell Composition of Stst.CTL_CKO, intAnno1", y = "Proportion") +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1, size = 9.6)) +
geom_point(aes(x=intAnno1, y = 100*prop, color= cnt),
position = position_dodge(0.75),
shape=16,alpha=0.75,size=2.15,
stroke=0.15, show.legend=F)
pcom1.intAnno1
stat1_intAnno2 <- test1.seur@meta.data[,c("cnt","rep","intAnno2")]
stat1_intAnno2.s <- stat1_intAnno2 %>%
group_by(cnt,rep,intAnno2) %>%
summarise(count=n()) %>%
mutate(prop= count/sum(count)) %>%
ungroup
## `summarise()` has grouped output by 'cnt', 'rep'. You can override using the
## `.groups` argument.
pcom1.intAnno2 <- stat1_intAnno2.s %>%
ggplot(aes(x = intAnno2, y = 100*prop, color=cnt)) +
geom_bar(stat="summary", fun="mean", position = position_dodge(0.75), width = 0.58, fill="white") +
theme_classic(base_size = 15) +
scale_color_manual(values = color.test1, name="") +
labs(title="Cell Composition of Stst.CTL_CKO, intAnno2", y = "Proportion") +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1, size = 9.6)) +
geom_point(aes(x=intAnno2, y = 100*prop, color= cnt),
position = position_dodge(0.75),
shape=16,alpha=0.75,size=2.15,
stroke=0.15, show.legend=F)
pcom1.intAnno2
glm.nb - anova.LRT
Sort.N <- c("Stst.CTL","Stst.CKO")
glm.nb_anova.1.intAnno1 <- list()
for(x1 in 1:2){
for(x2 in 1:2){
if(x2>x1){
stat1_intAnno1.s_N <- stat1_intAnno1.s %>% filter(cnt %in% c(Sort.N[x1],Sort.N[x2]))
total.N <- stat1_intAnno1.s_N %>%
group_by(cnt, rep) %>%
summarise(total=sum(count)) %>% ungroup() %>% as.data.frame()
rownames(total.N) <- paste0(as.character(total.N$cnt),
"_",
as.character(total.N$rep))
stat1_intAnno1.s_N$total <- total.N[paste0(stat1_intAnno1.s_N$cnt,
"_",
stat1_intAnno1.s_N$rep),"total"]
glm.nb_anova.1.intAnno1[[paste0(Sort.N[x1],"vs",Sort.N[x2])]] <- list()
for(AA in levels(stat1_intAnno1.s_N$intAnno1)){
glm.nb_anova.1.intAnno1[[paste0(Sort.N[x1],"vs",Sort.N[x2])]][[AA]] <- tryCatch({
aaa <- stat1_intAnno1.s_N %>% filter(intAnno1==AA)
aaa.glmnb <- MASS::glm.nb(count ~ cnt + offset(log(total)), data = aaa, maxit=1000)
aaa.anova <- anova(aaa.glmnb, test = "LRT")
aaa.anova$'Pr(>Chi)'[2]
}, error=function(e){
tryCatch({
aaa <- stat1_intAnno1.s_N %>% filter(intAnno1==AA)
aaa.glmnb <- MASS::glm.nb(count ~ cnt + offset(log(total)), data = aaa, maxit=100)
aaa.anova <- anova(aaa.glmnb, test = "LRT")
aaa.anova$'Pr(>Chi)'[2]
}, error=function(e){
tryCatch({
aaa <- stat1_intAnno1.s_N %>% filter(intAnno1==AA)
aaa.glmnb <- MASS::glm.nb(count ~ cnt + offset(log(total)), data = aaa, maxit=10)
aaa.anova <- anova(aaa.glmnb, test = "LRT")
aaa.anova$'Pr(>Chi)'[2]
}, error = function(e){
NA
})
})
})
}
glm.nb_anova.1.intAnno1[[paste0(Sort.N[x1],"vs",Sort.N[x2])]] <- unlist(glm.nb_anova.1.intAnno1[[paste0(Sort.N[x1],"vs",Sort.N[x2])]])
}
}
}
glm.nb_anova.1.intAnno1_df <- t(data.frame(Reduce(function(x,y){rbind(x,y)},glm.nb_anova.1.intAnno1)))
rownames(glm.nb_anova.1.intAnno1_df) <- names(glm.nb_anova.1.intAnno1)
colnames(glm.nb_anova.1.intAnno1_df) <- gsub("X","C",colnames(glm.nb_anova.1.intAnno1_df))
glm.nb_anova.1.intAnno1_df
## EMN1 EMN2 EMN3 EMN4 EMN5 IMN1
## Stst.CTLvsStst.CKO 0.437101 0.02226962 0.8008005 0.03400082 0.8765206 0.3068567
## IMN2 IMN3 IMN4 IN1 IN2
## Stst.CTLvsStst.CKO 0.002090672 0.2576225 0.07350664 2.335817e-05 0.4113903
## IN3 IPAN1 IPAN2 IPAN3 IPAN4
## Stst.CTLvsStst.CKO 0.06840589 0.7111369 0.8422819 0.4009551 0.888078
#options (scipen = 999)
round(glm.nb_anova.1.intAnno1_df,6)
## EMN1 EMN2 EMN3 EMN4 EMN5 IMN1
## Stst.CTLvsStst.CKO 0.437101 0.02227 0.800801 0.034001 0.876521 0.306857
## IMN2 IMN3 IMN4 IN1 IN2 IN3 IPAN1
## Stst.CTLvsStst.CKO 0.002091 0.257623 0.073507 2.3e-05 0.41139 0.068406 0.711137
## IPAN2 IPAN3 IPAN4
## Stst.CTLvsStst.CKO 0.842282 0.400955 0.888078
Sort.N <- c("Stst.CTL","Stst.CKO")
glm.nb_anova.1.intAnno2 <- list()
for(x1 in 1:2){
for(x2 in 1:2){
if(x2>x1){
stat1_intAnno2.s_N <- stat1_intAnno2.s %>% filter(cnt %in% c(Sort.N[x1],Sort.N[x2]))
total.N <- stat1_intAnno2.s_N %>%
group_by(cnt, rep) %>%
summarise(total=sum(count)) %>% ungroup() %>% as.data.frame()
rownames(total.N) <- paste0(as.character(total.N$cnt),
"_",
as.character(total.N$rep))
stat1_intAnno2.s_N$total <- total.N[paste0(stat1_intAnno2.s_N$cnt,
"_",
stat1_intAnno2.s_N$rep),"total"]
glm.nb_anova.1.intAnno2[[paste0(Sort.N[x1],"vs",Sort.N[x2])]] <- list()
for(AA in levels(stat1_intAnno2.s_N$intAnno2)){
glm.nb_anova.1.intAnno2[[paste0(Sort.N[x1],"vs",Sort.N[x2])]][[AA]] <- tryCatch({
aaa <- stat1_intAnno2.s_N %>% filter(intAnno2==AA)
aaa.glmnb <- MASS::glm.nb(count ~ cnt + offset(log(total)), data = aaa, maxit=1000)
aaa.anova <- anova(aaa.glmnb, test = "LRT")
aaa.anova$'Pr(>Chi)'[2]
}, error=function(e){
tryCatch({
aaa <- stat1_intAnno2.s_N %>% filter(intAnno2==AA)
aaa.glmnb <- MASS::glm.nb(count ~ cnt + offset(log(total)), data = aaa, maxit=100)
aaa.anova <- anova(aaa.glmnb, test = "LRT")
aaa.anova$'Pr(>Chi)'[2]
}, error=function(e){
tryCatch({
aaa <- stat1_intAnno2.s_N %>% filter(intAnno2==AA)
aaa.glmnb <- MASS::glm.nb(count ~ cnt + offset(log(total)), data = aaa, maxit=10)
aaa.anova <- anova(aaa.glmnb, test = "LRT")
aaa.anova$'Pr(>Chi)'[2]
}, error = function(e){
NA
})
})
})
}
glm.nb_anova.1.intAnno2[[paste0(Sort.N[x1],"vs",Sort.N[x2])]] <- unlist(glm.nb_anova.1.intAnno2[[paste0(Sort.N[x1],"vs",Sort.N[x2])]])
}
}
}
glm.nb_anova.1.intAnno2_df <- t(data.frame(Reduce(function(x,y){rbind(x,y)},glm.nb_anova.1.intAnno2)))
rownames(glm.nb_anova.1.intAnno2_df) <- names(glm.nb_anova.1.intAnno2)
colnames(glm.nb_anova.1.intAnno2_df) <- gsub("X","C",colnames(glm.nb_anova.1.intAnno2_df))
glm.nb_anova.1.intAnno2_df
## EMN1 EMN2 EMN3 EMN4 EMN5 IMN1
## Stst.CTLvsStst.CKO 0.437101 0.02226962 0.8008005 0.03400082 0.8765206 0.3068567
## IMN2 IMN3 IMN4 IN1 IN2
## Stst.CTLvsStst.CKO 0.002090672 0.2576225 0.07350664 2.335817e-05 0.4113903
## IN3 IPAN1.1 IPAN1.2 IPAN2.1 IPAN2.2 IPAN3
## Stst.CTLvsStst.CKO 0.06840589 0.2661605 0.8564025 0.5426043 0.5770597 0.4009551
## IPAN4
## Stst.CTLvsStst.CKO 0.888078
#options (scipen = 999)
round(glm.nb_anova.1.intAnno2_df,5)
## EMN1 EMN2 EMN3 EMN4 EMN5 IMN1 IMN2 IMN3
## Stst.CTLvsStst.CKO 0.4371 0.02227 0.8008 0.034 0.87652 0.30686 0.00209 0.25762
## IMN4 IN1 IN2 IN3 IPAN1.1 IPAN1.2 IPAN2.1
## Stst.CTLvsStst.CKO 0.07351 2e-05 0.41139 0.06841 0.26616 0.8564 0.5426
## IPAN2.2 IPAN3 IPAN4
## Stst.CTLvsStst.CKO 0.57706 0.40096 0.88808
cowplot::plot_grid(
pheatmap::pheatmap(table(cnt=test2.seur$sample,
clusters=test2.seur$intAnno1)[c(4:7,1:3),],
main = "Cell Count",
gaps_row = c(4),
gaps_col = c(5,9,12),
cluster_rows = F, cluster_cols = F, display_numbers = T, number_format = "%.0f", legend = F, silent = T)$gtable,
pheatmap::pheatmap(100*rowRatio(table(cnt=test2.seur$sample,
clusters=test2.seur$intAnno1)[c(4:7,1:3),] ),
main = "Cell Ratio",
gaps_row = c(4),
gaps_col = c(5,9,12),
cluster_rows = F, cluster_cols = F, display_numbers = T, number_format = "%.1f", legend = F, silent =T)$gtable,
ncol = 1)
cowplot::plot_grid(
pheatmap::pheatmap(table(cnt=test2.seur$sample,
clusters=test2.seur$intAnno2)[c(4:7,1:3),],
main = "Cell Count",
gaps_row = c(4),
gaps_col = c(5,9,12),
cluster_rows = F, cluster_cols = F, display_numbers = T, number_format = "%.0f", legend = F, silent = T)$gtable,
pheatmap::pheatmap(100*rowRatio(table(cnt=test2.seur$sample,
clusters=test2.seur$intAnno2)[c(4:7,1:3),] ),
main = "Cell Ratio",
gaps_row = c(4),
gaps_col = c(5,9,12),
cluster_rows = F, cluster_cols = F, display_numbers = T, number_format = "%.1f", legend = F, silent =T)$gtable,
ncol = 1)
stat2_intAnno1 <- test2.seur@meta.data[,c("cnt","rep","intAnno1")]
stat2_intAnno1.s <- stat2_intAnno1 %>%
group_by(cnt,rep,intAnno1) %>%
summarise(count=n()) %>%
mutate(prop= count/sum(count)) %>%
ungroup
## `summarise()` has grouped output by 'cnt', 'rep'. You can override using the
## `.groups` argument.
pcom2.intAnno1 <- stat2_intAnno1.s %>%
ggplot(aes(x = intAnno1, y = 100*prop, color=cnt)) +
geom_bar(stat="summary", fun="mean", position = position_dodge(0.75), width = 0.58, fill="white") +
theme_classic(base_size = 15) +
scale_color_manual(values = color.test2, name="") +
labs(title="Cell Composition of CR7d.CTL_CKO, intAnno1", y = "Proportion") +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1, size = 9.6)) +
geom_point(aes(x=intAnno1, y = 100*prop, color= cnt),
position = position_dodge(0.75),
shape=16,alpha=0.75,size=2.15,
stroke=0.15, show.legend=F)
pcom2.intAnno1
stat2_intAnno2 <- test2.seur@meta.data[,c("cnt","rep","intAnno2")]
stat2_intAnno2.s <- stat2_intAnno2 %>%
group_by(cnt,rep,intAnno2) %>%
summarise(count=n()) %>%
mutate(prop= count/sum(count)) %>%
ungroup
## `summarise()` has grouped output by 'cnt', 'rep'. You can override using the
## `.groups` argument.
pcom2.intAnno2 <- stat2_intAnno2.s %>%
ggplot(aes(x = intAnno2, y = 100*prop, color=cnt)) +
geom_bar(stat="summary", fun="mean", position = position_dodge(0.75), width = 0.58, fill="white") +
theme_classic(base_size = 15) +
scale_color_manual(values = color.test2, name="") +
labs(title="Cell Composition of CR7d.CTL_CKO, intAnno2", y = "Proportion") +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1, size = 9.6)) +
geom_point(aes(x=intAnno2, y = 100*prop, color= cnt),
position = position_dodge(0.75),
shape=16,alpha=0.75,size=2.15,
stroke=0.15, show.legend=F)
pcom2.intAnno2
glm.nb - anova.LRT
Sort.N <- c("CR7d.CTL","CR7d.CKO")
glm.nb_anova.2.intAnno1 <- list()
for(x1 in 1:2){
for(x2 in 1:2){
if(x2>x1){
stat2_intAnno1.s_N <- stat2_intAnno1.s %>% filter(cnt %in% c(Sort.N[x1],Sort.N[x2]))
total.N <- stat2_intAnno1.s_N %>%
group_by(cnt, rep) %>%
summarise(total=sum(count)) %>% ungroup() %>% as.data.frame()
rownames(total.N) <- paste0(as.character(total.N$cnt),
"_",
as.character(total.N$rep))
stat2_intAnno1.s_N$total <- total.N[paste0(stat2_intAnno1.s_N$cnt,
"_",
stat2_intAnno1.s_N$rep),"total"]
glm.nb_anova.2.intAnno1[[paste0(Sort.N[x1],"vs",Sort.N[x2])]] <- list()
for(AA in levels(stat2_intAnno1.s_N$intAnno1)){
glm.nb_anova.2.intAnno1[[paste0(Sort.N[x1],"vs",Sort.N[x2])]][[AA]] <- tryCatch({
aaa <- stat2_intAnno1.s_N %>% filter(intAnno1==AA)
aaa.glmnb <- MASS::glm.nb(count ~ cnt + offset(log(total)), data = aaa, maxit=1000)
aaa.anova <- anova(aaa.glmnb, test = "LRT")
aaa.anova$'Pr(>Chi)'[2]
}, error=function(e){
tryCatch({
aaa <- stat2_intAnno1.s_N %>% filter(intAnno1==AA)
aaa.glmnb <- MASS::glm.nb(count ~ cnt + offset(log(total)), data = aaa, maxit=100)
aaa.anova <- anova(aaa.glmnb, test = "LRT")
aaa.anova$'Pr(>Chi)'[2]
}, error=function(e){
tryCatch({
aaa <- stat2_intAnno1.s_N %>% filter(intAnno1==AA)
aaa.glmnb <- MASS::glm.nb(count ~ cnt + offset(log(total)), data = aaa, maxit=10)
aaa.anova <- anova(aaa.glmnb, test = "LRT")
aaa.anova$'Pr(>Chi)'[2]
}, error = function(e){
NA
})
})
})
}
glm.nb_anova.2.intAnno1[[paste0(Sort.N[x1],"vs",Sort.N[x2])]] <- unlist(glm.nb_anova.2.intAnno1[[paste0(Sort.N[x1],"vs",Sort.N[x2])]])
}
}
}
glm.nb_anova.2.intAnno1_df <- t(data.frame(Reduce(function(x,y){rbind(x,y)},glm.nb_anova.2.intAnno1)))
rownames(glm.nb_anova.2.intAnno1_df) <- names(glm.nb_anova.2.intAnno1)
colnames(glm.nb_anova.2.intAnno1_df) <- gsub("X","C",colnames(glm.nb_anova.2.intAnno1_df))
glm.nb_anova.2.intAnno1_df
## EMN1 EMN2 EMN3 EMN4 EMN5 IMN1
## CR7d.CTLvsCR7d.CKO 0.3964044 0.4170971 0.5165039 0.6549232 0.469626 0.2860499
## IMN2 IMN3 IMN4 IN1 IN2
## CR7d.CTLvsCR7d.CKO 0.01985923 0.5701387 0.9166744 0.001347195 0.7900928
## IN3 IPAN1 IPAN2 IPAN3 IPAN4
## CR7d.CTLvsCR7d.CKO 0.4605807 0.1626049 0.6988405 0.09326914 0.285167
#options (scipen = 999)
round(glm.nb_anova.2.intAnno1_df,4)
## EMN1 EMN2 EMN3 EMN4 EMN5 IMN1 IMN2 IMN3
## CR7d.CTLvsCR7d.CKO 0.3964 0.4171 0.5165 0.6549 0.4696 0.286 0.0199 0.5701
## IMN4 IN1 IN2 IN3 IPAN1 IPAN2 IPAN3 IPAN4
## CR7d.CTLvsCR7d.CKO 0.9167 0.0013 0.7901 0.4606 0.1626 0.6988 0.0933 0.2852
Sort.N <- c("CR7d.CTL","CR7d.CKO")
glm.nb_anova.2.intAnno2 <- list()
for(x1 in 1:2){
for(x2 in 1:2){
if(x2>x1){
stat2_intAnno2.s_N <- stat2_intAnno2.s %>% filter(cnt %in% c(Sort.N[x1],Sort.N[x2]))
total.N <- stat2_intAnno2.s_N %>%
group_by(cnt, rep) %>%
summarise(total=sum(count)) %>% ungroup() %>% as.data.frame()
rownames(total.N) <- paste0(as.character(total.N$cnt),
"_",
as.character(total.N$rep))
stat2_intAnno2.s_N$total <- total.N[paste0(stat2_intAnno2.s_N$cnt,
"_",
stat2_intAnno2.s_N$rep),"total"]
glm.nb_anova.2.intAnno2[[paste0(Sort.N[x1],"vs",Sort.N[x2])]] <- list()
for(AA in levels(stat2_intAnno2.s_N$intAnno2)){
glm.nb_anova.2.intAnno2[[paste0(Sort.N[x1],"vs",Sort.N[x2])]][[AA]] <- tryCatch({
aaa <- stat2_intAnno2.s_N %>% filter(intAnno2==AA)
aaa.glmnb <- MASS::glm.nb(count ~ cnt + offset(log(total)), data = aaa, maxit=1000)
aaa.anova <- anova(aaa.glmnb, test = "LRT")
aaa.anova$'Pr(>Chi)'[2]
}, error=function(e){
tryCatch({
aaa <- stat2_intAnno2.s_N %>% filter(intAnno2==AA)
aaa.glmnb <- MASS::glm.nb(count ~ cnt + offset(log(total)), data = aaa, maxit=100)
aaa.anova <- anova(aaa.glmnb, test = "LRT")
aaa.anova$'Pr(>Chi)'[2]
}, error=function(e){
tryCatch({
aaa <- stat2_intAnno2.s_N %>% filter(intAnno2==AA)
aaa.glmnb <- MASS::glm.nb(count ~ cnt + offset(log(total)), data = aaa, maxit=10)
aaa.anova <- anova(aaa.glmnb, test = "LRT")
aaa.anova$'Pr(>Chi)'[2]
}, error = function(e){
NA
})
})
})
}
glm.nb_anova.2.intAnno2[[paste0(Sort.N[x1],"vs",Sort.N[x2])]] <- unlist(glm.nb_anova.2.intAnno2[[paste0(Sort.N[x1],"vs",Sort.N[x2])]])
}
}
}
glm.nb_anova.2.intAnno2_df <- t(data.frame(Reduce(function(x,y){rbind(x,y)},glm.nb_anova.2.intAnno2)))
rownames(glm.nb_anova.2.intAnno2_df) <- names(glm.nb_anova.2.intAnno2)
colnames(glm.nb_anova.2.intAnno2_df) <- gsub("X","C",colnames(glm.nb_anova.2.intAnno2_df))
glm.nb_anova.2.intAnno2_df
## EMN1 EMN2 EMN3 EMN4 EMN5 IMN1
## CR7d.CTLvsCR7d.CKO 0.3964044 0.4170971 0.5165039 0.6549232 0.469626 0.2860499
## IMN2 IMN3 IMN4 IN1 IN2
## CR7d.CTLvsCR7d.CKO 0.01985923 0.5701387 0.9166744 0.001347195 0.7900928
## IN3 IPAN1.1 IPAN1.2 IPAN2.1 IPAN2.2 IPAN3
## CR7d.CTLvsCR7d.CKO 0.4605807 0.1757535 0.3257703 0.9831953 0.5659296 0.09326914
## IPAN4
## CR7d.CTLvsCR7d.CKO 0.285167
#options (scipen = 999)
round(glm.nb_anova.2.intAnno2_df,4)
## EMN1 EMN2 EMN3 EMN4 EMN5 IMN1 IMN2 IMN3
## CR7d.CTLvsCR7d.CKO 0.3964 0.4171 0.5165 0.6549 0.4696 0.286 0.0199 0.5701
## IMN4 IN1 IN2 IN3 IPAN1.1 IPAN1.2 IPAN2.1 IPAN2.2
## CR7d.CTLvsCR7d.CKO 0.9167 0.0013 0.7901 0.4606 0.1758 0.3258 0.9832 0.5659
## IPAN3 IPAN4
## CR7d.CTLvsCR7d.CKO 0.0933 0.2852
# find markers for every cluster compared to all remaining cells
Idents(GEX.seur) <- "intAnno1"
GEX.seur <- PrepSCTFindMarkers(GEX.seur, assay = "SCT")
#GEX.markers.pre <- FindAllMarkers(GEX.seur, only.pos = TRUE, min.pct = 0.01,
# assay = "SCT",
# test.use = "MAST",
# logfc.threshold = 0.1)
GEX.markers.pre <- read.table("Baf53cre_CR.markers.SCT_intAnno1.202405.csv", header = TRUE, sep = ",")
#GEX.markers.pre %>% group_by(cluster) %>% top_n(n = 8, wt = avg_log2FC)
GEX.markers.pre$cluster <- factor(as.character(GEX.markers.pre$cluster),
levels = levels(GEX.seur$intAnno1))
markers.pre_t60 <- (GEX.markers.pre %>% group_by(cluster) %>%
filter(pct.1>0.05 & gene %in% grep("Rps|Rpl|mt-|Hsp",gene,invert = T,value = T)) %>%
top_n(n = 60, wt = avg_log2FC) %>%
ungroup() %>%
arrange(desc(avg_log2FC*pct.1),gene) %>%
distinct(gene, .keep_all = TRUE) %>%
arrange(cluster,p_val_adj))$gene
markers.pre_t120 <- (GEX.markers.pre %>% group_by(cluster) %>%
filter(pct.1>0.05 & gene %in% grep("Rps|Rpl|mt-|Hsp",gene,invert = T,value = T)) %>%
top_n(n = 120, wt = avg_log2FC) %>%
ungroup() %>%
arrange(desc(avg_log2FC*pct.1),gene) %>%
distinct(gene, .keep_all = TRUE) %>%
arrange(cluster,p_val_adj))$gene
ttt = 468
ttt/60
## [1] 7.8
ttt/64
## [1] 7.3125
ttt/65
## [1] 7.2
pp.t60 <- list()
for(i in 1:8){
pp.t60[[i]] <- DotPlot(GEX.seur, features = rev(markers.pre_t60[(60*i-59):(60*i)])) + coord_flip() +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1, size = 9.6))
}
pp.t60
## [[1]]
##
## [[2]]
##
## [[3]]
##
## [[4]]
##
## [[5]]
##
## [[6]]
##
## [[7]]
##
## [[8]]
ttt = 839
ttt/60
## [1] 13.98333
ttt/64
## [1] 13.10938
ttt/65
## [1] 12.90769
pp.t120 <- list()
for(i in 1:14){
pp.t120[[i]] <- DotPlot(GEX.seur, features = rev(markers.pre_t120[(60*i-59):(60*i)])) + coord_flip() +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1, size = 9.6))
}
pp.t120
## [[1]]
##
## [[2]]
##
## [[3]]
##
## [[4]]
##
## [[5]]
##
## [[6]]
##
## [[7]]
##
## [[8]]
##
## [[9]]
##
## [[10]]
##
## [[11]]
##
## [[12]]
##
## [[13]]
##
## [[14]]
head(GEX.seur@meta.data)
## orig.ident nCount_RNA nFeature_RNA percent.mt percent.rb
## AAACCCAAGAATACAC-1_2 Stst.CTL_CKO 2282 1290 0.1314636 0.1314636
## AAACCCAAGCAATAGT-1_2 Stst.CTL_CKO 1535 1011 0.1302932 0.5863192
## AAACCCAAGGTGAGCT-1_2 Stst.CTL_CKO 1049 737 0.0000000 0.4766444
## AAACCCACAACGAGGT-1_2 Stst.CTL_CKO 1591 1046 0.6913891 0.4399749
## AAACCCACAAGAGTAT-1_2 Stst.CTL_CKO 2106 1250 0.2849003 0.2849003
## AAACCCACAATCAAGA-1_2 Stst.CTL_CKO 1451 953 0.1378360 0.1378360
## FB.info S.Score G2M.Score Phase cnt rep
## AAACCCAAGAATACAC-1_2 CTL.3 -0.00110357 0.012058898 G2M Stst.CTL rep3
## AAACCCAAGCAATAGT-1_2 CKO.3 0.02907907 -0.015558699 S Stst.CKO rep3
## AAACCCAAGGTGAGCT-1_2 CTL.1 -0.01100811 -0.008015087 G1 Stst.CTL rep1
## AAACCCACAACGAGGT-1_2 CTL.1 -0.02143685 0.005086498 G2M Stst.CTL rep1
## AAACCCACAAGAGTAT-1_2 CKO.3 0.02965845 -0.009057774 S Stst.CKO rep3
## AAACCCACAATCAAGA-1_2 CTL.2 -0.01216686 0.004143546 G2M Stst.CTL rep2
## sample tissue nCount_SCT nFeature_SCT seurat_clusters
## AAACCCAAGAATACAC-1_2 Stst.CTL3 Ileum 1865 1287 8
## AAACCCAAGCAATAGT-1_2 Stst.CKO3 Ileum 1525 1010 2
## AAACCCAAGGTGAGCT-1_2 Stst.CTL1 Ileum 1170 735 0
## AAACCCACAACGAGGT-1_2 Stst.CTL1 Ileum 1571 1045 3
## AAACCCACAAGAGTAT-1_2 Stst.CKO3 Ileum 1837 1248 16
## AAACCCACAATCAAGA-1_2 Stst.CTL2 Ileum 1452 953 3
## sort_clusters intAnno1 intAnno2
## AAACCCAAGAATACAC-1_2 8 EMN1 EMN1
## AAACCCAAGCAATAGT-1_2 2 IMN1 IMN1
## AAACCCAAGGTGAGCT-1_2 0 EMN1 EMN1
## AAACCCACAACGAGGT-1_2 3 EMN1 EMN1
## AAACCCACAAGAGTAT-1_2 16 IMN3 IMN3
## AAACCCACAATCAAGA-1_2 3 EMN1 EMN1
Idents(GEX.seur) <- "cnt"
#GEX.seur <- PrepSCTFindMarkers(GEX.seur, assay = "SCT")
neur.clusters <- grep("EMN|IMN|IN|IPAN",levels(GEX.seur$intAnno2), value = T)
neur.clusters
## [1] "EMN1" "EMN2" "EMN3" "EMN4" "EMN5" "IMN1" "IMN2"
## [8] "IMN3" "IMN4" "IN1" "IN2" "IN3" "IPAN1.1" "IPAN1.2"
## [15] "IPAN2.1" "IPAN2.2" "IPAN3" "IPAN4"
#
all_new <- neur.clusters
all_new
## [1] "EMN1" "EMN2" "EMN3" "EMN4" "EMN5" "IMN1" "IMN2"
## [8] "IMN3" "IMN4" "IN1" "IN2" "IN3" "IPAN1.1" "IPAN1.2"
## [15] "IPAN2.1" "IPAN2.2" "IPAN3" "IPAN4"
#
list_new <- list()
list_new[["All"]] <- all_new
list_new[["EMNs"]] <- grep("EMN",all_new,value = T)
list_new[["IMNs"]] <- grep("IMN",all_new,value = T)
list_new[["IPAN1"]] <- grep("IPAN1",all_new,value = T)
list_new[["IPAN2"]] <- grep("IPAN2",all_new,value = T)
#list_new[["INs"]] <- grep("IN",all_new,value = T)
#list_new[["IPANs"]] <- grep("IPAN",all_new,value = T)
for(NN in all_new){
list_new[[NN]] <- NN
}
names_new <- names(list_new)
list_new
## $All
## [1] "EMN1" "EMN2" "EMN3" "EMN4" "EMN5" "IMN1" "IMN2"
## [8] "IMN3" "IMN4" "IN1" "IN2" "IN3" "IPAN1.1" "IPAN1.2"
## [15] "IPAN2.1" "IPAN2.2" "IPAN3" "IPAN4"
##
## $EMNs
## [1] "EMN1" "EMN2" "EMN3" "EMN4" "EMN5"
##
## $IMNs
## [1] "IMN1" "IMN2" "IMN3" "IMN4"
##
## $IPAN1
## [1] "IPAN1.1" "IPAN1.2"
##
## $IPAN2
## [1] "IPAN2.1" "IPAN2.2"
##
## $EMN1
## [1] "EMN1"
##
## $EMN2
## [1] "EMN2"
##
## $EMN3
## [1] "EMN3"
##
## $EMN4
## [1] "EMN4"
##
## $EMN5
## [1] "EMN5"
##
## $IMN1
## [1] "IMN1"
##
## $IMN2
## [1] "IMN2"
##
## $IMN3
## [1] "IMN3"
##
## $IMN4
## [1] "IMN4"
##
## $IN1
## [1] "IN1"
##
## $IN2
## [1] "IN2"
##
## $IN3
## [1] "IN3"
##
## $IPAN1.1
## [1] "IPAN1.1"
##
## $IPAN1.2
## [1] "IPAN1.2"
##
## $IPAN2.1
## [1] "IPAN2.1"
##
## $IPAN2.2
## [1] "IPAN2.2"
##
## $IPAN3
## [1] "IPAN3"
##
## $IPAN4
## [1] "IPAN4"
names_new
## [1] "All" "EMNs" "IMNs" "IPAN1" "IPAN2" "EMN1" "EMN2"
## [8] "EMN3" "EMN4" "EMN5" "IMN1" "IMN2" "IMN3" "IMN4"
## [15] "IN1" "IN2" "IN3" "IPAN1.1" "IPAN1.2" "IPAN2.1" "IPAN2.2"
## [22] "IPAN3" "IPAN4"
#test1.seur <- subset(GEX.seur, subset= cnt %in% c("Stst.CTL","Stst.CKO"))
test1.seur
## An object of class Seurat
## 47749 features across 11409 samples within 3 assays
## Active assay: SCT (21008 features, 0 variable features)
## 2 other assays present: RNA, integrated
## 3 dimensional reductions calculated: pca, tsne, umap
#head(test1.seur@active.ident)
Idents(test1.seur) <- "cnt"
#test1.DEGs_new
# save DEGs' table
df_test1.DEGs_new <- data.frame()
for(nn in names_new){
df_test1.DEGs_new <- rbind(df_test1.DEGs_new, data.frame(test1.DEGs_new[[nn]],intAnno=nn))
}
#write.csv(df_test1.DEGs_new, "Baf53cre_CR.DEGs.Stst_CTLvsCKO.intAnno2.csv")
#write.csv(df_test1.DEGs_new %>% filter(intAnno %in% setdiff(names_new,paste0("IPAN",c(1.1,1.2,2.1,2.2)))), "Baf53cre_CR.DEGs.Stst_CTLvsCKO.intAnno1.csv")
cut0: padj 0.05 log2FC 0.1
cut1: padj 0.05 log2FC 0.3
cut2: padj 0.01 log2FC > log2(1.5) (|FC|>1.5)
# cut1
cut.padj = 0.05
cut.log2FC = 0.1
cut.pct1 = 0.1
stat_test1o.DEGs_new <- df_test1.DEGs_new %>% filter(p_val_adj < cut.padj &
abs(avg_log2FC) > cut.log2FC &
pct.1 > cut.pct1) %>%
group_by(intAnno,cluster) %>%
summarise(up.DEGs = n()) %>% as.data.frame()
stat_test1o.DEGs_new
## intAnno cluster up.DEGs
## 1 All Stst.CTL 15
## 2 All Stst.CKO 17
## 3 EMN1 Stst.CTL 7
## 4 EMN1 Stst.CKO 4
## 5 EMN2 Stst.CTL 1
## 6 EMN3 Stst.CTL 1
## 7 EMN4 Stst.CTL 1
## 8 EMN4 Stst.CKO 1
## 9 EMN5 Stst.CTL 1
## 10 EMNs Stst.CTL 11
## 11 EMNs Stst.CKO 8
## 12 IMN1 Stst.CTL 4
## 13 IMN1 Stst.CKO 2
## 14 IMN2 Stst.CTL 1
## 15 IMN3 Stst.CTL 1
## 16 IMNs Stst.CTL 6
## 17 IMNs Stst.CKO 5
## 18 IN1 Stst.CTL 2
## 19 IN2 Stst.CTL 1
## 20 IPAN1 Stst.CTL 4
## 21 IPAN1 Stst.CKO 2
## 22 IPAN1.1 Stst.CTL 1
## 23 IPAN1.1 Stst.CKO 1
## 24 IPAN1.2 Stst.CTL 2
## 25 IPAN2 Stst.CTL 3
## 26 IPAN2 Stst.CKO 2
## 27 IPAN2.1 Stst.CTL 2
## 28 IPAN2.1 Stst.CKO 2
## 29 IPAN2.2 Stst.CTL 1
## 30 IPAN4 Stst.CTL 1
## 31 IPAN4 Stst.CKO 1
df_test1.DEGs_new %>% filter(p_val_adj < cut.padj &
abs(avg_log2FC) > cut.log2FC &
pct.1 > cut.pct1) %>%
group_by(intAnno,cluster) %>%
summarise(up.DEGs = n()) %>% as.data.frame() %>%
ggplot(aes(x=intAnno, y=up.DEGs, color = cluster)) +
geom_bar(stat="summary", fun="mean", position = position_dodge(0.75), width = 0.58, fill="white") +
geom_text(aes(label = up.DEGs),vjust=-0.21, show.legend = F, position = position_dodge(0.75) ) +
theme_classic(base_size = 15) +
scale_color_manual(values = color.test1, name="") +
labs(title=paste0("up.DEGs stat, pct.1>",cut.pct1,", padj<",cut.padj,", |log2FC|>",cut.log2FC), y = "Count") +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1, size = 9.6),
title =element_text(size=12, face='bold'))
# cut1
cut.padj = 0.05
cut.log2FC = 0.3
cut.pct1 = 0.1
stat_test1a.DEGs_new <- df_test1.DEGs_new %>% filter(p_val_adj < cut.padj &
abs(avg_log2FC) > cut.log2FC &
pct.1 > cut.pct1) %>%
group_by(intAnno,cluster) %>%
summarise(up.DEGs = n()) %>% as.data.frame()
stat_test1a.DEGs_new
## intAnno cluster up.DEGs
## 1 All Stst.CTL 1
## 2 All Stst.CKO 2
## 3 EMN1 Stst.CTL 4
## 4 EMN1 Stst.CKO 3
## 5 EMN2 Stst.CTL 1
## 6 EMN3 Stst.CTL 1
## 7 EMN4 Stst.CTL 1
## 8 EMN4 Stst.CKO 1
## 9 EMN5 Stst.CTL 1
## 10 EMNs Stst.CTL 2
## 11 EMNs Stst.CKO 3
## 12 IMN1 Stst.CTL 1
## 13 IMN1 Stst.CKO 1
## 14 IMN2 Stst.CTL 1
## 15 IMN3 Stst.CTL 1
## 16 IMNs Stst.CTL 1
## 17 IMNs Stst.CKO 2
## 18 IN1 Stst.CTL 2
## 19 IN2 Stst.CTL 1
## 20 IPAN1 Stst.CTL 4
## 21 IPAN1 Stst.CKO 2
## 22 IPAN1.1 Stst.CTL 1
## 23 IPAN1.1 Stst.CKO 1
## 24 IPAN1.2 Stst.CTL 2
## 25 IPAN2 Stst.CTL 3
## 26 IPAN2 Stst.CKO 2
## 27 IPAN2.1 Stst.CTL 2
## 28 IPAN2.1 Stst.CKO 2
## 29 IPAN2.2 Stst.CTL 1
## 30 IPAN4 Stst.CTL 1
## 31 IPAN4 Stst.CKO 1
df_test1.DEGs_new %>% filter(p_val_adj < cut.padj &
abs(avg_log2FC) > cut.log2FC &
pct.1 > cut.pct1) %>%
group_by(intAnno,cluster) %>%
summarise(up.DEGs = n()) %>% as.data.frame() %>%
ggplot(aes(x=intAnno, y=up.DEGs, color = cluster)) +
geom_bar(stat="summary", fun="mean", position = position_dodge(0.75), width = 0.58, fill="white") +
geom_text(aes(label = up.DEGs),vjust=-0.21, show.legend = F, position = position_dodge(0.75) ) +
theme_classic(base_size = 15) +
scale_color_manual(values = color.test1, name="") +
labs(title=paste0("up.DEGs stat, pct.1>",cut.pct1,", padj<",cut.padj,", |log2FC|>",cut.log2FC), y = "Count") +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1, size = 9.6),
title =element_text(size=12, face='bold'))
# cut2
cut.padj = 0.01
cut.log2FC = log2(1.5)
cut.pct1 = 0.1
stat_test1b.DEGs_new <- df_test1.DEGs_new %>% filter(p_val_adj < cut.padj &
abs(avg_log2FC) > cut.log2FC &
pct.1 > cut.pct1) %>%
group_by(intAnno,cluster) %>%
summarise(up.DEGs = n()) %>% as.data.frame()
stat_test1b.DEGs_new
## intAnno cluster up.DEGs
## 1 All Stst.CTL 1
## 2 EMN1 Stst.CTL 1
## 3 EMN2 Stst.CTL 1
## 4 EMN3 Stst.CTL 1
## 5 EMN4 Stst.CTL 1
## 6 EMN4 Stst.CKO 1
## 7 EMN5 Stst.CTL 1
## 8 EMNs Stst.CTL 1
## 9 IMN1 Stst.CTL 1
## 10 IMN2 Stst.CTL 1
## 11 IMN3 Stst.CTL 1
## 12 IMNs Stst.CTL 1
## 13 IN1 Stst.CTL 2
## 14 IPAN1 Stst.CTL 2
## 15 IPAN1 Stst.CKO 1
## 16 IPAN1.1 Stst.CTL 1
## 17 IPAN1.1 Stst.CKO 1
## 18 IPAN1.2 Stst.CTL 1
## 19 IPAN2 Stst.CTL 2
## 20 IPAN2.1 Stst.CTL 2
## 21 IPAN4 Stst.CTL 1
df_test1.DEGs_new %>% filter(p_val_adj < cut.padj &
abs(avg_log2FC) > cut.log2FC &
pct.1 > cut.pct1) %>%
group_by(intAnno,cluster) %>%
summarise(up.DEGs = n()) %>% as.data.frame() %>%
ggplot(aes(x=intAnno, y=up.DEGs, color = cluster)) +
geom_bar(stat="summary", fun="mean", position = position_dodge(0.75), width = 0.58, fill="white") +
geom_text(aes(label = up.DEGs),vjust=-0.21, show.legend = F, position = position_dodge(0.75) ) +
theme_classic(base_size = 15) +
scale_color_manual(values = color.test1, name="") +
labs(title=paste0("up.DEGs stat, pct.1>",cut.pct1,", padj<",cut.padj,", |FC|>",2^cut.log2FC), y = "Count") +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1, size = 9.6),
title =element_text(size=12, face='bold'))
pp_test1.DEGs.new <- lapply(list_new, function(x){NA})
#
test1.DEGs_new.combine <- test1.DEGs_new
test1.DEGs_new.combine <- lapply(test1.DEGs_new.combine, function(x){
x[x$cluster=="Stst.CTL","avg_log2FC"] <- -x[x$cluster=="Stst.CTL","avg_log2FC"]
x
})
pp_test1.DEGs.new$All <- test1.DEGs_new.combine$All %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="All Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$All
pp_test1.DEGs.new$EMNs <- test1.DEGs_new.combine$EMNs %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="EMNs Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$EMNs
pp_test1.DEGs.new$EMN1 <- test1.DEGs_new.combine$EMN1 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="EMN1 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$EMN1
pp_test1.DEGs.new$EMN2 <- test1.DEGs_new.combine$EMN2 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="EMN2 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$EMN2
pp_test1.DEGs.new$EMN3 <- test1.DEGs_new.combine$EMN3 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="EMN3 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$EMN3
pp_test1.DEGs.new$EMN4 <- test1.DEGs_new.combine$EMN4 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="EMN4 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$EMN4
pp_test1.DEGs.new$EMN5 <- test1.DEGs_new.combine$EMN5 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="EMN5 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$EMN5
pp_test1.DEGs.new$IMNs <- test1.DEGs_new.combine$IMNs %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="IMNs Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$IMNs
pp_test1.DEGs.new$IMN1 <- test1.DEGs_new.combine$IMN1 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="IMN1 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$IMN1
pp_test1.DEGs.new$IMN2 <- test1.DEGs_new.combine$IMN2 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="IMN2 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$IMN2
pp_test1.DEGs.new$IMN3 <- test1.DEGs_new.combine$IMN3 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="IMN3 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$IMN3
pp_test1.DEGs.new$IMN4 <- test1.DEGs_new.combine$IMN4 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="IMN4 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$IMN4
pp_test1.DEGs.new$IN1 <- test1.DEGs_new.combine$IN1 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="IN1 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$IN1
pp_test1.DEGs.new$IN2 <- test1.DEGs_new.combine$IN2 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="IN2 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$IN2
pp_test1.DEGs.new$IN3 <- test1.DEGs_new.combine$IN3 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="IN3 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$IN3
pp_test1.DEGs.new$IPAN1 <- test1.DEGs_new.combine$IPAN1 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="IPAN1 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$IPAN1
pp_test1.DEGs.new$IPAN1.1 <- test1.DEGs_new.combine$IPAN1.1 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="IPAN1.1 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$IPAN1.1
pp_test1.DEGs.new$IPAN1.2 <- test1.DEGs_new.combine$IPAN1.2 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="IPAN1.2 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$IPAN1.2
pp_test1.DEGs.new$IPAN2 <- test1.DEGs_new.combine$IPAN2 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="IPAN2 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$IPAN2
pp_test1.DEGs.new$IPAN2.1 <- test1.DEGs_new.combine$IPAN2.1 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="IPAN2.1 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$IPAN2.1
pp_test1.DEGs.new$IPAN2.2 <- test1.DEGs_new.combine$IPAN2.2 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="IPAN2.2 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$IPAN2.2
pp_test1.DEGs.new$IPAN3 <- test1.DEGs_new.combine$IPAN3 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="IPAN3 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$IPAN3
pp_test1.DEGs.new$IPAN4 <- test1.DEGs_new.combine$IPAN4 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"Stst.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"Stst.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("Stst.CTL"=as.vector(color.test1[1]),
"Stst.CKO"=as.vector(color.test1[2]),
"None"="grey")) +
theme_classic() + labs(title="IPAN4 Stst.CTL vs Stst.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test1.DEGs.new$IPAN4
#test2.seur <- subset(GEX.seur, subset= cnt %in% c("CR7d.CTL","CR7d.CKO"))
test2.seur
## An object of class Seurat
## 47749 features across 12899 samples within 3 assays
## Active assay: SCT (21008 features, 0 variable features)
## 2 other assays present: RNA, integrated
## 3 dimensional reductions calculated: pca, tsne, umap
#head(test2.seur@active.ident)
Idents(test2.seur) <- "cnt"
#test2.DEGs_new
# save DEGs' table
df_test2.DEGs_new <- data.frame()
for(nn in names_new){
df_test2.DEGs_new <- rbind(df_test2.DEGs_new, data.frame(test2.DEGs_new[[nn]],intAnno=nn))
}
#write.csv(df_test2.DEGs_new, "Baf53cre_CR.DEGs.CR7d_CTLvsCKO.intAnno2.csv")
#write.csv(df_test2.DEGs_new %>% filter(intAnno %in% setdiff(names_new,paste0("IPAN",c(1.1,1.2,2.1,2.2)))), "Baf53cre_CR.DEGs.CR7d_CTLvsCKO.intAnno1.csv")
cut0: padj 0.05 log2FC 0.1
cut1: padj 0.05 log2FC 0.3
cut2: padj 0.01 log2FC > log2(1.5) (|FC|>1.5)
# cut1
cut.padj = 0.05
cut.log2FC = 0.1
cut.pct1 = 0.1
stat_test2o.DEGs_new <- df_test2.DEGs_new %>% filter(p_val_adj < cut.padj &
abs(avg_log2FC) > cut.log2FC &
pct.1 > cut.pct1) %>%
group_by(intAnno,cluster) %>%
summarise(up.DEGs = n()) %>% as.data.frame()
stat_test2o.DEGs_new
## intAnno cluster up.DEGs
## 1 All CR7d.CTL 4
## 2 All CR7d.CKO 12
## 3 EMN1 CR7d.CTL 3
## 4 EMN1 CR7d.CKO 3
## 5 EMN2 CR7d.CTL 2
## 6 EMN3 CR7d.CTL 2
## 7 EMN4 CR7d.CTL 2
## 8 EMN4 CR7d.CKO 2
## 9 EMN5 CR7d.CTL 2
## 10 EMNs CR7d.CTL 3
## 11 EMNs CR7d.CKO 6
## 12 IMN1 CR7d.CTL 3
## 13 IMN1 CR7d.CKO 4
## 14 IMN2 CR7d.CTL 2
## 15 IMN3 CR7d.CTL 1
## 16 IMNs CR7d.CTL 4
## 17 IMNs CR7d.CKO 4
## 18 IN1 CR7d.CTL 1
## 19 IN1 CR7d.CKO 1
## 20 IPAN1 CR7d.CTL 1
## 21 IPAN1 CR7d.CKO 1
## 22 IPAN1.1 CR7d.CTL 1
## 23 IPAN1.2 CR7d.CTL 1
## 24 IPAN1.2 CR7d.CKO 1
## 25 IPAN2 CR7d.CTL 1
## 26 IPAN2 CR7d.CKO 3
## 27 IPAN2.1 CR7d.CTL 1
## 28 IPAN2.1 CR7d.CKO 1
## 29 IPAN2.2 CR7d.CTL 1
## 30 IPAN2.2 CR7d.CKO 1
df_test2.DEGs_new %>% filter(p_val_adj < cut.padj &
abs(avg_log2FC) > cut.log2FC &
pct.1 > cut.pct1) %>%
group_by(intAnno,cluster) %>%
summarise(up.DEGs = n()) %>% as.data.frame() %>%
ggplot(aes(x=intAnno, y=up.DEGs, color = cluster)) +
geom_bar(stat="summary", fun="mean", position = position_dodge(0.75), width = 0.58, fill="white") +
geom_text(aes(label = up.DEGs),vjust=-0.21, show.legend = F, position = position_dodge(0.75) ) +
theme_classic(base_size = 15) +
scale_color_manual(values = color.test2, name="") +
labs(title=paste0("up.DEGs stat, pct.1>",cut.pct1,", padj<",cut.padj,", |log2FC|>",cut.log2FC), y = "Count") +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1, size = 9.6),
title =element_text(size=12, face='bold'))
# cut1
cut.padj = 0.05
cut.log2FC = 0.3
cut.pct1 = 0.1
stat_test2a.DEGs_new <- df_test2.DEGs_new %>% filter(p_val_adj < cut.padj &
abs(avg_log2FC) > cut.log2FC &
pct.1 > cut.pct1) %>%
group_by(intAnno,cluster) %>%
summarise(up.DEGs = n()) %>% as.data.frame()
stat_test2a.DEGs_new
## intAnno cluster up.DEGs
## 1 All CR7d.CTL 1
## 2 All CR7d.CKO 2
## 3 EMN1 CR7d.CTL 3
## 4 EMN1 CR7d.CKO 2
## 5 EMN2 CR7d.CTL 1
## 6 EMN3 CR7d.CTL 2
## 7 EMN4 CR7d.CTL 2
## 8 EMN4 CR7d.CKO 2
## 9 EMN5 CR7d.CTL 2
## 10 EMNs CR7d.CTL 3
## 11 EMNs CR7d.CKO 2
## 12 IMN1 CR7d.CTL 1
## 13 IMN1 CR7d.CKO 4
## 14 IMN2 CR7d.CTL 2
## 15 IMN3 CR7d.CTL 1
## 16 IMNs CR7d.CTL 2
## 17 IMNs CR7d.CKO 4
## 18 IN1 CR7d.CTL 1
## 19 IN1 CR7d.CKO 1
## 20 IPAN1 CR7d.CTL 1
## 21 IPAN1 CR7d.CKO 1
## 22 IPAN1.1 CR7d.CTL 1
## 23 IPAN1.2 CR7d.CTL 1
## 24 IPAN1.2 CR7d.CKO 1
## 25 IPAN2 CR7d.CTL 1
## 26 IPAN2 CR7d.CKO 3
## 27 IPAN2.1 CR7d.CTL 1
## 28 IPAN2.1 CR7d.CKO 1
## 29 IPAN2.2 CR7d.CTL 1
## 30 IPAN2.2 CR7d.CKO 1
df_test2.DEGs_new %>% filter(p_val_adj < cut.padj &
abs(avg_log2FC) > cut.log2FC &
pct.1 > cut.pct1) %>%
group_by(intAnno,cluster) %>%
summarise(up.DEGs = n()) %>% as.data.frame() %>%
ggplot(aes(x=intAnno, y=up.DEGs, color = cluster)) +
geom_bar(stat="summary", fun="mean", position = position_dodge(0.75), width = 0.58, fill="white") +
geom_text(aes(label = up.DEGs),vjust=-0.21, show.legend = F, position = position_dodge(0.75) ) +
theme_classic(base_size = 15) +
scale_color_manual(values = color.test2, name="") +
labs(title=paste0("up.DEGs stat, pct.1>",cut.pct1,", padj<",cut.padj,", |log2FC|>",cut.log2FC), y = "Count") +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1, size = 9.6),
title =element_text(size=12, face='bold'))
# cut2
cut.padj = 0.01
cut.log2FC = log2(1.5)
cut.pct1 = 0.1
stat_test2b.DEGs_new <- df_test2.DEGs_new %>% filter(p_val_adj < cut.padj &
abs(avg_log2FC) > cut.log2FC &
pct.1 > cut.pct1) %>%
group_by(intAnno,cluster) %>%
summarise(up.DEGs = n()) %>% as.data.frame()
stat_test2b.DEGs_new
## intAnno cluster up.DEGs
## 1 All CR7d.CTL 1
## 2 EMN1 CR7d.CTL 1
## 3 EMN2 CR7d.CTL 1
## 4 EMN3 CR7d.CTL 1
## 5 EMN4 CR7d.CTL 1
## 6 EMN5 CR7d.CTL 2
## 7 EMNs CR7d.CTL 1
## 8 IMN1 CR7d.CTL 1
## 9 IMN2 CR7d.CTL 1
## 10 IMNs CR7d.CTL 1
## 11 IN1 CR7d.CTL 1
## 12 IN1 CR7d.CKO 1
## 13 IPAN1 CR7d.CTL 1
## 14 IPAN1 CR7d.CKO 1
## 15 IPAN1.1 CR7d.CTL 1
## 16 IPAN1.2 CR7d.CTL 1
## 17 IPAN1.2 CR7d.CKO 1
## 18 IPAN2 CR7d.CTL 1
## 19 IPAN2 CR7d.CKO 2
## 20 IPAN2.1 CR7d.CTL 1
## 21 IPAN2.2 CR7d.CTL 1
## 22 IPAN2.2 CR7d.CKO 1
df_test2.DEGs_new %>% filter(p_val_adj < cut.padj &
abs(avg_log2FC) > cut.log2FC &
pct.1 > cut.pct1) %>%
group_by(intAnno,cluster) %>%
summarise(up.DEGs = n()) %>% as.data.frame() %>%
ggplot(aes(x=intAnno, y=up.DEGs, color = cluster)) +
geom_bar(stat="summary", fun="mean", position = position_dodge(0.75), width = 0.58, fill="white") +
geom_text(aes(label = up.DEGs),vjust=-0.21, show.legend = F, position = position_dodge(0.75) ) +
theme_classic(base_size = 15) +
scale_color_manual(values = color.test2, name="") +
labs(title=paste0("up.DEGs stat, pct.1>",cut.pct1,", padj<",cut.padj,", |FC|>",2^cut.log2FC), y = "Count") +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1, size = 9.6),
title =element_text(size=12, face='bold'))
pp_test2.DEGs.new <- lapply(list_new, function(x){NA})
#
test2.DEGs_new.combine <- test2.DEGs_new
test2.DEGs_new.combine <- lapply(test2.DEGs_new.combine, function(x){
x[x$cluster=="CR7d.CTL","avg_log2FC"] <- -x[x$cluster=="CR7d.CTL","avg_log2FC"]
x
})
pp_test2.DEGs.new$All <- test2.DEGs_new.combine$All %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="All CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$All
pp_test2.DEGs.new$EMNs <- test2.DEGs_new.combine$EMNs %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="EMNs CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$EMNs
pp_test2.DEGs.new$EMN1 <- test2.DEGs_new.combine$EMN1 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="EMN1 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$EMN1
pp_test2.DEGs.new$EMN2 <- test2.DEGs_new.combine$EMN2 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="EMN2 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$EMN2
pp_test2.DEGs.new$EMN3 <- test2.DEGs_new.combine$EMN3 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="EMN3 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$EMN3
pp_test2.DEGs.new$EMN4 <- test2.DEGs_new.combine$EMN4 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="EMN4 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$EMN4
pp_test2.DEGs.new$EMN5 <- test2.DEGs_new.combine$EMN5 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="EMN5 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$EMN5
pp_test2.DEGs.new$IMNs <- test2.DEGs_new.combine$IMNs %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="IMNs CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$IMNs
pp_test2.DEGs.new$IMN1 <- test2.DEGs_new.combine$IMN1 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="IMN1 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$IMN1
pp_test2.DEGs.new$IMN2 <- test2.DEGs_new.combine$IMN2 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="IMN2 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$IMN2
pp_test2.DEGs.new$IMN3 <- test2.DEGs_new.combine$IMN3 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="IMN3 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$IMN3
pp_test2.DEGs.new$IMN4 <- test2.DEGs_new.combine$IMN4 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="IMN4 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$IMN4
pp_test2.DEGs.new$IN1 <- test2.DEGs_new.combine$IN1 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="IN1 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$IN1
pp_test2.DEGs.new$IN2 <- test2.DEGs_new.combine$IN2 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="IN2 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$IN2
pp_test2.DEGs.new$IN3 <- test2.DEGs_new.combine$IN3 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="IN3 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$IN3
pp_test2.DEGs.new$IPAN1 <- test2.DEGs_new.combine$IPAN1 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="IPAN1 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$IPAN1
pp_test2.DEGs.new$IPAN1.1 <- test2.DEGs_new.combine$IPAN1.1 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="IPAN1.1 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$IPAN1.1
pp_test2.DEGs.new$IPAN1.2 <- test2.DEGs_new.combine$IPAN1.2 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="IPAN1.2 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$IPAN1.2
pp_test2.DEGs.new$IPAN2 <- test2.DEGs_new.combine$IPAN2 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="IPAN2 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$IPAN2
pp_test2.DEGs.new$IPAN2.1 <- test2.DEGs_new.combine$IPAN2.1 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="IPAN2.1 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$IPAN2.1
pp_test2.DEGs.new$IPAN2.2 <- test2.DEGs_new.combine$IPAN2.2 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="IPAN2.2 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$IPAN2.2
pp_test2.DEGs.new$IPAN3 <- test2.DEGs_new.combine$IPAN3 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="IPAN3 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$IPAN3
pp_test2.DEGs.new$IPAN4 <- test2.DEGs_new.combine$IPAN4 %>%
mutate(label=ifelse(((p_val_adj < 1e-3 & avg_log2FC<0) | (p_val_adj < 1e-3 & avg_log2FC>0) | (p_val_adj < 0.05 & abs(avg_log2FC) > 0.1)) &
#!grepl("^Gm|^Hsp|^Rps|^Rpl|Rik$|Xist|Tsix|Uty|Ddx3y|Eif2s3y|Kdm5d",gene),gene,""),
!grepl("^Rps|^Rpl",gene),gene,""),
sig.=ifelse(p_val_adj<0.05 & avg_log2FC<0,"CR7d.CTL",ifelse(p_val_adj<0.05 & avg_log2FC>0,"CR7d.CKO","None")),
padj=ifelse(p_val_adj<1e-60,p_val_adj+1e-60,p_val_adj)) %>%
ggplot(aes(y= -log10(padj), x=avg_log2FC, size=pct.1, label=label,fill=sig.)) +
geom_point(shape=21, alpha=0.65) +
geom_text_repel(size=2.5, max.overlaps = 500) +
scale_fill_manual(values = c("CR7d.CTL"=as.vector(color.test2[1]),
"CR7d.CKO"=as.vector(color.test2[2]),
"None"="grey")) +
theme_classic() + labs(title="IPAN4 CR7d.CTL vs CR7d.CKO") +
guides(size = guide_legend(order = 2), fill=guide_legend(order = 1, reverse = F)) +
geom_vline(xintercept = c(-0.1,0.1), linetype="dotted")
pp_test2.DEGs.new$IPAN4
markers.old.s <- list(EMN=c("Chat","Bnc2",#"Tox","Ptprt",
"Gfra2","Oprk1",#"Adamtsl1",
"Fbxw15","Fbxw24",#"Chrna7",
"Satb1","Cntnap5b",
"Gabrb1","Nxph1","Lama2","Lrrc7",
"Ryr3",#"Eda",
"Tac1",
#"Kctd8","Ntrk2",
"Penk",
"Fut9","Nfatc1","Egfr","Ahr"#,#"Mgll",
#"Chrm3"
),
IMN=c("Nos1","Kcnab1",
"Gfra1","Etv1",
#"Man1a","Airn",
"Adcy2","Cmah","Creb5","Vip","Pde1a",
"Ebf1"#,"Gpc5"
),
IN=c("Npas3","Synpr","St18","Gal",
"Neurod6",
#"Kcnk13",
"Moxd1","Sctr",
"Piezo1","Sst",#"Adamts9",
"Kcnn2"),
IPAN=c("Calb2","Calcb","Nmu","Adgrg6",#"Pcdh10",
"Ngfr","Galr1","Il7",#"Aff2",
#"Gpr149",
"Cdh6","Cdh8",
"Clstn2",#"Ano2","Ntrk3",
"Cpne4",#"Vwc2l",
"Cdh9","Scgn",
#"Vcan",
"Cck","Piezo2","Kcnh7",
#"Rerg",
"Bmpr1b","Skap1","Ntng1",
"Tafa2","Nxph2"))
pn.intAnno1.test0a <- DotPlot(GEX.seur, features = as.vector(unlist(markers.old.s)), group.by = "intAnno1",
cols = c("midnightblue","darkorange1")) +
#coord_flip() +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) +
scale_y_discrete(limits=rev) +
labs(title="All intAnno1")
pn.intAnno1.test0a
pn.intAnno2.test0a <- DotPlot(GEX.seur, features = as.vector(unlist(markers.old.s)), group.by = "intAnno2",
cols = c("midnightblue","darkorange1")) +
#coord_flip() +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) +
scale_y_discrete(limits=rev) +
labs(title="All intAnno2")
pn.intAnno2.test0a
pn.intAnno1.test0b <- DotPlot(GEX.seur, features = as.vector(unlist(markers.old.s)), group.by = "intAnno1",
cols = c("midnightblue","darkorange1")) +
coord_flip() +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.4,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) + scale_y_discrete(limits=rev) +
labs(title="All intAnno1")
pn.intAnno1.test0b
pn.intAnno2.test0b <- DotPlot(GEX.seur, features = as.vector(unlist(markers.old.s)), group.by = "intAnno2",
cols = c("midnightblue","darkorange1")) +
coord_flip() +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.4,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) + scale_y_discrete(limits=rev) +
labs(title="All intAnno2")
pn.intAnno2.test0b
pn.intAnno1.test1a <- DotPlot(test1.seur, features = as.vector(unlist(markers.old.s)), group.by = "intAnno1",
cols = c("midnightblue","darkorange1")) +
#coord_flip() +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) +
scale_y_discrete(limits=rev) +
labs(title="Stst intAnno1")
pn.intAnno1.test1a
pn.intAnno2.test1a <- DotPlot(test1.seur, features = as.vector(unlist(markers.old.s)), group.by = "intAnno2",
cols = c("midnightblue","darkorange1")) +
#coord_flip() +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) +
scale_y_discrete(limits=rev) +
labs(title="Stst intAnno2")
pn.intAnno2.test1a
pn.intAnno1.test1b <- DotPlot(test1.seur, features = as.vector(unlist(markers.old.s)), group.by = "intAnno1",
cols = c("midnightblue","darkorange1")) +
coord_flip() +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.4,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) + scale_y_discrete(limits=rev) +
labs(title="Stst intAnno1")
pn.intAnno1.test1b
pn.intAnno2.test1b <- DotPlot(test1.seur, features = as.vector(unlist(markers.old.s)), group.by = "intAnno2",
cols = c("midnightblue","darkorange1")) +
coord_flip() +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.4,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) + scale_y_discrete(limits=rev) +
labs(title="Stst intAnno2")
pn.intAnno2.test1b
pn.intAnno1.test2a <- DotPlot(test2.seur, features = as.vector(unlist(markers.old.s)), group.by = "intAnno1",
cols = c("midnightblue","darkorange1")) +
#coord_flip() +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) +
scale_y_discrete(limits=rev) +
labs(title="CR7d intAnno1")
pn.intAnno1.test2a
pn.intAnno2.test2a <- DotPlot(test2.seur, features = as.vector(unlist(markers.old.s)), group.by = "intAnno2",
cols = c("midnightblue","darkorange1")) +
#coord_flip() +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) +
scale_y_discrete(limits=rev) +
labs(title="CR7d intAnno2")
pn.intAnno2.test2a
pn.intAnno1.test2b <- DotPlot(test2.seur, features = as.vector(unlist(markers.old.s)), group.by = "intAnno1",
cols = c("midnightblue","darkorange1")) +
coord_flip() +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.4,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) + scale_y_discrete(limits=rev) +
labs(title="CR7d intAnno1")
pn.intAnno1.test2b
pn.intAnno2.test2b <- DotPlot(test2.seur, features = as.vector(unlist(markers.old.s)), group.by = "intAnno2",
cols = c("midnightblue","darkorange1")) +
coord_flip() +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.4,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) + scale_y_discrete(limits=rev) +
labs(title="CR7d intAnno2")
pn.intAnno2.test2b
ref.seur <- readRDS("../../20230704_10x_SZJ/analysis_ref/GSE149524.P21.integration_Anno.s.rds")
ref.seur
## An object of class Seurat
## 37583 features across 4419 samples within 3 assays
## Active assay: SCT (16365 features, 0 variable features)
## 2 other assays present: RNA, integrated
## 3 dimensional reductions calculated: pca, tsne, umap
DimPlot(ref.seur, reduction = "umap", label = T, group.by = "Anno1", cols = color.ref) +
DimPlot(ref.seur, reduction = "umap", label = T, group.by = "Anno2")
pn.ref.a <- DotPlot(ref.seur, features = as.vector(unlist(markers.old.s)), group.by = "Anno2",
cols = c("midnightblue","darkorange1")) +
#coord_flip() +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) +
scale_y_discrete(limits=rev) +
labs(title="NatNeur2021 P21")
pn.ref.a
pn.ref.b <- DotPlot(ref.seur, features = as.vector(unlist(markers.old.s)), group.by = "Anno2",
cols = c("midnightblue","darkorange1")) +
coord_flip() +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.4,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) + scale_y_discrete(limits=rev) +
labs(title="NatNeur2021 P21")
pn.ref.b
markers.new.s <- list(EMN=c("Chat","Bnc2",#"Tox","Ptprt",
"Gfra2","Oprk1",#"Adamtsl1",
"Fbxw15","Fbxw24",#"Chrna7",
"Satb1","Itga6","Cntnap5b",
"Chgb","Nxph1",
"Lama2","Efnb2","Itgb8",
"Lrrc7",
"Ryr3",#"Eda",
"Tac1",
#"Kctd8","Ntrk2",
"Penk",
"Fut9","Nfatc1","Egfr","Ahr"#"Mgll",
#"Chrm3"
),
IMN=c("Nos1","Kcnab1",
"Gfra1","Etv1",
#"Man1a","Airn",
"Adcy2","Cmah","Col25a1",
"Mid1","Creb5","Vip","Pde1a",
"Ebf1",#,"Gpc5"
"Ppara","Pcdh11x",
"Adcy8","Grp"
),
IN=c("Npas3","Synpr","St18","Gal",
"Cdh10","Neurod6",
"Kcnk13",
"Moxd1","Sctr",
"Piezo1","Vipr2","Sst",#"Adamts9",
"Kcnn2"
),
IPAN=c("Calb2","Adcy1",
"Nmu","Adgrg6",#"Pcdh10",
"Ngfr","Il7",
"Itgb6","Calcb","Galr1",
#"Aff2",
#"Gpr149",
"Met",
"Cpne4","Cdh6","Cdh8",
"Clstn2",#"Ano2","Ntrk3",
#"Vwc2l",
"Car10","Scgn","Glp2r","Cck",
"Cdh9",
#"Vcan",
"Dcc",
"Gabrb1",
"Piezo2","Kcnh7",
#"Rerg",
"Bmpr1b","Ntng1","Skap1",
"Tafa2","Nxph2"))
pm.intAnno1.test0a <- DotPlot(GEX.seur, features = as.vector(unlist(markers.new.s)), group.by = "intAnno1",
cols = c("midnightblue","darkorange1")) +
#coord_flip() +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) +
scale_y_discrete(limits=rev) +
labs(title="All intAnno1")
pm.intAnno1.test0a
pm.intAnno2.test0a <- DotPlot(GEX.seur, features = as.vector(unlist(markers.new.s)), group.by = "intAnno2",
cols = c("midnightblue","darkorange1")) +
#coord_flip() +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) +
scale_y_discrete(limits=rev) +
labs(title="All intAnno2")
pm.intAnno2.test0a
pm.intAnno1.test0b <- DotPlot(GEX.seur, features = as.vector(unlist(markers.new.s)), group.by = "intAnno1",
cols = c("midnightblue","darkorange1")) +
coord_flip() +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.4,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) + scale_y_discrete(limits=rev) +
labs(title="All intAnno1")
pm.intAnno1.test0b
pm.intAnno2.test0b <- DotPlot(GEX.seur, features = as.vector(unlist(markers.new.s)), group.by = "intAnno2",
cols = c("midnightblue","darkorange1")) +
coord_flip() +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.4,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) + scale_y_discrete(limits=rev) +
labs(title="All intAnno2")
pm.intAnno2.test0b
pm.intAnno1.test1a <- DotPlot(test1.seur, features = as.vector(unlist(markers.new.s)), group.by = "intAnno1",
cols = c("midnightblue","darkorange1")) +
#coord_flip() +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) +
scale_y_discrete(limits=rev) +
labs(title="Stst intAnno1")
pm.intAnno1.test1a
pm.intAnno2.test1a <- DotPlot(test1.seur, features = as.vector(unlist(markers.new.s)), group.by = "intAnno2",
cols = c("midnightblue","darkorange1")) +
#coord_flip() +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) +
scale_y_discrete(limits=rev) +
labs(title="Stst intAnno2")
pm.intAnno2.test1a
pm.intAnno1.test1b <- DotPlot(test1.seur, features = as.vector(unlist(markers.new.s)), group.by = "intAnno1",
cols = c("midnightblue","darkorange1")) +
coord_flip() +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.4,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) + scale_y_discrete(limits=rev) +
labs(title="Stst intAnno1")
pm.intAnno1.test1b
pm.intAnno2.test1b <- DotPlot(test1.seur, features = as.vector(unlist(markers.new.s)), group.by = "intAnno2",
cols = c("midnightblue","darkorange1")) +
coord_flip() +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.4,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) + scale_y_discrete(limits=rev) +
labs(title="Stst intAnno2")
pm.intAnno2.test1b
pm.intAnno1.test2a <- DotPlot(test2.seur, features = as.vector(unlist(markers.new.s)), group.by = "intAnno1",
cols = c("midnightblue","darkorange1")) +
#coord_flip() +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) +
scale_y_discrete(limits=rev) +
labs(title="CR7d intAnno1")
pm.intAnno1.test2a
pm.intAnno2.test2a <- DotPlot(test2.seur, features = as.vector(unlist(markers.new.s)), group.by = "intAnno2",
cols = c("midnightblue","darkorange1")) +
#coord_flip() +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) +
scale_y_discrete(limits=rev) +
labs(title="CR7d intAnno2")
pm.intAnno2.test2a
pm.intAnno1.test2b <- DotPlot(test2.seur, features = as.vector(unlist(markers.new.s)), group.by = "intAnno1",
cols = c("midnightblue","darkorange1")) +
coord_flip() +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.4,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) + scale_y_discrete(limits=rev) +
labs(title="CR7d intAnno1")
pm.intAnno1.test2b
pm.intAnno2.test2b <- DotPlot(test2.seur, features = as.vector(unlist(markers.new.s)), group.by = "intAnno2",
cols = c("midnightblue","darkorange1")) +
coord_flip() +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.4,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) + scale_y_discrete(limits=rev) +
labs(title="CR7d intAnno2")
pm.intAnno2.test2b
pm.ref.a <- DotPlot(ref.seur, features = as.vector(unlist(markers.new.s)), group.by = "Anno2",
cols = c("midnightblue","darkorange1")) +
#coord_flip() +
theme(axis.text.x = element_text(angle = 60, hjust = 1, vjust = 1,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) +
scale_y_discrete(limits=rev) +
labs(title="NatNeur2021 P21")
pm.ref.a
pm.ref.b <- DotPlot(ref.seur, features = as.vector(unlist(markers.new.s)), group.by = "Anno2",
cols = c("midnightblue","darkorange1")) +
coord_flip() +
theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.4,size = 9.15))+
scale_color_gradientn(colours = material.heat(100)) +
#scale_x_discrete(limits=rev) + scale_y_discrete(limits=rev) +
labs(title="NatNeur2021 P21")
pm.ref.b
use top120-built top markers
load("I:/Shared_win/projects/20231220_10x_SZJ/analysis_plus/top.markerlist.RData")
data.frame(row.names = "topmarkers",lapply(top.list,length))
## EMN1 EMN2 EMN3 EMN4 EMN5 IMN1 IMN2 IMN3 IMN4 IN1 IN2 IN3 IPAN1.1
## topmarkers 96 32 37 34 70 71 33 49 53 60 56 63 70
## IPAN1.2 IPAN2.1 IPAN2.2 IPAN3 IPAN4
## topmarkers 51 43 52 52 50
for(nn in neur.clusters){
GEX.seur <- add_geneset_score(obj = GEX.seur,
Assay = "SCT",
geneset = top.list[[nn]],
setname = paste0("score.",nn))
}
## Summarizing data
## Summarizing data
## Summarizing data
## Summarizing data
## Summarizing data
## Summarizing data
## Summarizing data
## Summarizing data
## Summarizing data
## Summarizing data
## Summarizing data
## Summarizing data
## Summarizing data
## Summarizing data
## Summarizing data
## Summarizing data
## Summarizing data
## Summarizing data
mapal <- colorRampPalette(RColorBrewer::brewer.pal(4,"Spectral"))(120)
FeaturePlot(GEX.seur, features = paste0("score.",neur.clusters), ncol = 5) &
scale_color_gradientn(colors = rev(mapal))
mapal <- colorRampPalette(RColorBrewer::brewer.pal(4,"Spectral"))(120)
FeaturePlot(GEX.seur, features = paste0("score.",neur.clusters), ncol = 5, raster = T, pt.size = 3.5) &
scale_color_gradientn(colors = rev(mapal))
check.m <- "Nmu"
#
VlnPlot(test1.seur, features = check.m, group.by = "intAnno1", assay = "RNA") + NoLegend() + labs(x="Stst") +
scale_fill_manual(values = color.A1)
VlnPlot(test1.seur, features = check.m, group.by = "intAnno1", assay = "RNA", split.by = "cnt") + NoLegend() + labs(x="Stst") +
scale_fill_manual(values = color.test1)
#
VlnPlot(test2.seur, features = check.m, group.by = "intAnno1", assay = "RNA") + NoLegend() + labs(x="CR7d") +
scale_fill_manual(values = color.A1)
VlnPlot(test2.seur, features = check.m, group.by = "intAnno1", assay = "RNA", split.by = "cnt") + NoLegend() + labs(x="CR7d") +
scale_fill_manual(values = color.test2)
##
# use coord_cartesian(ylim()) instead of ylim()
VlnPlot(test1.seur, features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("All - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,4.5)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(3.5),
size=3.5
)
#
VlnPlot(subset(test1.seur, subset=intAnno1 %in% c("IPAN1")), features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,4.5)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(3.5),
size=3.5
)
VlnPlot(subset(test1.seur, subset=intAnno2 %in% c("IPAN1.1")), features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1.1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,4.5)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(3.5),
size=3.5
)
VlnPlot(subset(test1.seur, subset=intAnno2 %in% c("IPAN1.2")), features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1.2 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,4.5)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(3.5),
size=3.5
)
##
# use coord_cartesian(ylim()) instead of ylim()
VlnPlot(test2.seur, features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("All - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,4.5)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(3.5),
size=3.5
)
#
VlnPlot(subset(test2.seur, subset=intAnno1 %in% c("IPAN1")), features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,4.5)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(3.5),
size=3.5
)
VlnPlot(subset(test2.seur, subset=intAnno2 %in% c("IPAN1.1")), features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1.1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,4.5)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(3.5),
size=3.5
)
VlnPlot(subset(test2.seur, subset=intAnno2 %in% c("IPAN1.2")), features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1.2 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,4.5)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(3.5),
size=3.5
)
check.m <- "Calcb"
#
cowplot::plot_grid(
VlnPlot(test1.seur, features = check.m, group.by = "intAnno1", assay = "RNA") + NoLegend() + labs(x="Stst") +
scale_fill_manual(values = color.A1),
VlnPlot(test2.seur, features = check.m, group.by = "intAnno1", assay = "RNA") + NoLegend() + labs(x="CR7d") +
scale_fill_manual(values = color.A1),
VlnPlot(test1.seur, features = check.m, group.by = "intAnno1", assay = "RNA", split.by = "cnt") + NoLegend() + labs(x="Stst") +
scale_fill_manual(values = color.test1),
VlnPlot(test2.seur, features = check.m, group.by = "intAnno1", assay = "RNA", split.by = "cnt") + NoLegend() + labs(x="CR7d") +
scale_fill_manual(values = color.test2),
#
ncol = 2)
##
# use coord_cartesian(ylim()) instead of ylim()
cowplot::plot_grid(
VlnPlot(test1.seur, features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("All - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,3.5)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(2.75),
size=3.5
),
#
VlnPlot(subset(test1.seur, subset=intAnno1 %in% c("IPAN1")), features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,3.5)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(2.75),
size=3.5
) ,
VlnPlot(subset(test1.seur, subset=intAnno2 %in% c("IPAN1.1")), features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1.1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,3.5)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(2.75),
size=3.5
),
VlnPlot(subset(test1.seur, subset=intAnno2 %in% c("IPAN1.2")), features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1.2 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,3.5)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(2.75),
size=3.5
),
ncol = 2)
##
# use coord_cartesian(ylim()) instead of ylim()
cowplot::plot_grid(
VlnPlot(test2.seur, features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("All - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,3.5)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(2.75),
size=3.5
),
#
VlnPlot(subset(test2.seur, subset=intAnno1 %in% c("IPAN1")), features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,3.5)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(2.75),
size=3.5
),
VlnPlot(subset(test2.seur, subset=intAnno2 %in% c("IPAN1.1")), features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1.1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,3.5)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(2.75),
size=3.5
),
VlnPlot(subset(test2.seur, subset=intAnno2 %in% c("IPAN1.2")), features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1.2 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,3.5)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(2.75),
size=3.5
),
ncol = 2)
##
# use coord_cartesian(ylim()) instead of ylim()
cowplot::plot_grid(
VlnPlot(test1.seur, features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "SCT") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("All - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,2)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(1.5),
size=3.5
),
#
VlnPlot(subset(test1.seur, subset=intAnno1 %in% c("IPAN1")), features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "SCT") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,2)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(1.5),
size=3.5
) ,
VlnPlot(subset(test1.seur, subset=intAnno2 %in% c("IPAN1.1")), features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "SCT") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1.1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,2)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(1.5),
size=3.5
),
VlnPlot(subset(test1.seur, subset=intAnno2 %in% c("IPAN1.2")), features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "SCT") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1.2 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,2)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(1.5),
size=3.5
),
ncol = 2)
##
# use coord_cartesian(ylim()) instead of ylim()
cowplot::plot_grid(
VlnPlot(test2.seur, features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "SCT") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("All - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,2)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(1.45),
size=3.5
),
#
VlnPlot(subset(test2.seur, subset=intAnno1 %in% c("IPAN1")), features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "SCT") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,2)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(1.45),
size=3.5
),
VlnPlot(subset(test2.seur, subset=intAnno2 %in% c("IPAN1.1")), features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "SCT") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1.1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,2)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(1.25),
size=3.5
),
VlnPlot(subset(test2.seur, subset=intAnno2 %in% c("IPAN1.2")), features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "SCT") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1.2 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,2)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(1.45),
size=3.5
),
ncol = 2)
check.m <- "Gal"
#
cowplot::plot_grid(
VlnPlot(test1.seur, features = check.m, group.by = "intAnno1", assay = "RNA") + NoLegend() + labs(x="Stst") +
scale_fill_manual(values = color.A1),
VlnPlot(test2.seur, features = check.m, group.by = "intAnno1", assay = "RNA") + NoLegend() + labs(x="CR7d") +
scale_fill_manual(values = color.A1),
VlnPlot(test1.seur, features = check.m, group.by = "intAnno1", assay = "RNA", split.by = "cnt") + NoLegend() + labs(x="Stst") +
scale_fill_manual(values = color.test1),
VlnPlot(test2.seur, features = check.m, group.by = "intAnno1", assay = "RNA", split.by = "cnt") + NoLegend() + labs(x="CR7d") +
scale_fill_manual(values = color.test2),
#
ncol = 2)
##
# use coord_cartesian(ylim()) instead of ylim()
cowplot::plot_grid(
VlnPlot(test1.seur, features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("All - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(4.5),
size=3.5
),
#
VlnPlot(subset(test1.seur, subset=intAnno1 %in% c("IN1")), features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IN1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(4.5),
size=3.5
) ,
VlnPlot(subset(test1.seur, subset=intAnno2 %in% c("IN2")), features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IN2 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(4.5),
size=3.5
),
VlnPlot(subset(test1.seur, subset=intAnno2 %in% c("IMN3")), features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IMN3 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(2.75),
size=3.5
),
ncol = 2)
##
# use coord_cartesian(ylim()) instead of ylim()
cowplot::plot_grid(
VlnPlot(test2.seur, features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("All - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(4.5),
size=3.5
),
#
VlnPlot(subset(test2.seur, subset=intAnno1 %in% c("IN1")), features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IN1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(4.5),
size=3.5
) ,
VlnPlot(subset(test2.seur, subset=intAnno2 %in% c("IN2")), features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IN2 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(4.5),
size=3.5
),
VlnPlot(subset(test2.seur, subset=intAnno2 %in% c("IMN3")), features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IMN3 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(2.75),
size=3.5
),
ncol = 2)
check.m <- "Nrg1"
#
cowplot::plot_grid(
VlnPlot(test1.seur, features = check.m, group.by = "intAnno1", assay = "RNA") + NoLegend() + labs(x="Stst") +
scale_fill_manual(values = color.A1),
VlnPlot(test2.seur, features = check.m, group.by = "intAnno1", assay = "RNA") + NoLegend() + labs(x="CR7d") +
scale_fill_manual(values = color.A1),
VlnPlot(test1.seur, features = check.m, group.by = "intAnno1", assay = "RNA", split.by = "cnt") + NoLegend() + labs(x="Stst") +
scale_fill_manual(values = color.test1),
VlnPlot(test2.seur, features = check.m, group.by = "intAnno1", assay = "RNA", split.by = "cnt") + NoLegend() + labs(x="CR7d") +
scale_fill_manual(values = color.test2),
#
ncol = 2)
##
# use coord_cartesian(ylim()) instead of ylim()
cowplot::plot_grid(
VlnPlot(test1.seur, features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("All - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(4.5),
size=3.5
),
#
VlnPlot(subset(test1.seur, subset=intAnno1 %in% c("IN1")), features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IN1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(4.5),
size=3.5
) ,
VlnPlot(subset(test1.seur, subset=intAnno1 %in% c("IN2")), features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IN2 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(4),
size=3.5
),
VlnPlot(subset(test1.seur, subset=intAnno1 %in% c("IN3")), features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IN3 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(4),
size=3.5
),
VlnPlot(subset(test1.seur, subset=intAnno1 %in% c("IPAN1")), features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(4.5),
size=3.5
),
VlnPlot(subset(test1.seur, subset=intAnno1 %in% c("IMN4")), features = check.m, group.by = "cnt", cols = color.test1, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IMN4 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("Stst.CTL","Stst.CKO")),
label.y = c(3.25),
size=3.5
),
ncol = 2)
##
# use coord_cartesian(ylim()) instead of ylim()
cowplot::plot_grid(
VlnPlot(test2.seur, features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("All - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(4.5),
size=3.5
),
#
VlnPlot(subset(test2.seur, subset=intAnno1 %in% c("IN1")), features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IN1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(4.5),
size=3.5
) ,
VlnPlot(subset(test2.seur, subset=intAnno1 %in% c("IN2")), features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IN2 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(4),
size=3.5
),
VlnPlot(subset(test2.seur, subset=intAnno1 %in% c("IN3")), features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IN3 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(4),
size=3.5
),
VlnPlot(subset(test2.seur, subset=intAnno1 %in% c("IPAN1")), features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IPAN1 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(4.5),
size=3.5
),
VlnPlot(subset(test2.seur, subset=intAnno1 %in% c("IMN4")), features = check.m, group.by = "cnt", cols = color.test2, pt.size = 0, assay = "RNA") + NoLegend() &
geom_jitter(alpha=0.25, shape=16, width = 0.3, size = 0.25) & labs(title = paste0("IMN4 - ",check.m), x="") &
geom_boxplot(outlier.size = 0, fill="white", width=0.2, size=0.1, alpha=0.55) &
stat_summary(fun=mean, geom="point", shape=18, size=3, color="black", alpha=0.55) & coord_cartesian(ylim = c(0,6)) &
ggpubr::stat_compare_means(aes(lable = ..p.signif..),
method = "wilcox.test",
comparisons = list(c("CR7d.CTL","CR7d.CKO")),
label.y = c(3.25),
size=3.5
),
ncol = 2)